Exact solution of the $2d$ dimer model: Corner free energy, correlation functions and combinatorics

In this work, some classical results of the pfaffian theory of the dimer model based on the work of Kasteleyn, Fisher and Temperley are introduced in a fermionic framework. Then we shall detail the bosonic formulation of the model {\it via} the so-called height mapping and the nature of boundary conditions is unravelled. The complete and detailed fermionic solution of the dimer model on the square lattice with an arbitrary number of monomers is presented, and finite size effect analysis is performed to study surface and corner effects, leading to the extrapolation of the central charge of the model. The solution allows for exact calculations of monomer and dimer correlation functions in the discrete level and the scaling behavior can be inferred in order to find the set of scaling dimensions and compare to the bosonic theory which predict particular features concerning corner behaviors. Finally, some combinatorial and numerical properties of partition functions with boundary monomers are discussed, proved and checked with enumeration algorithms.

Perfect matching of the square lattice, and (b) its "domino" representation. This combinatorial problem reduces to the calculation of the partition function (1) with t = 1.

Introduction
Following Onsager's solution of the 2d Ising model in the forties [113,80,81], the introduction of the Bethe ansatz [13] and the discovery of the machinery of transfer matrices, the field of exact solutions of lattice statistical physics models has exploded leading to the birth of a new domain of theoretical and mathematical physics known as exactly solved models [10]. The confluence of this new field with 2d conformal field theory (CFT) [11] discovered by Belavin, Polyakov and Zamolodchikov (see [35] for an extensive monograph) had a huge impact in theoretical physics, from high energy to condensed matter, leading to a whole new level of understanding of classical and quantum integrable systems [97].
Initially, the dimer model (see Fig. 1) has been introduced by physicists to describe absorption of diatomic molecules on a 2d substrate [49], yet it became quickly a general problem studied in various scientific communities. From the mathematical point of view, this problem known as perfect matching problem [105]-is a famous and active problem of combinatorics and graph theory [47] with a large spectrum of applications. The enumeration of so-called Kekulé structures of molecular graphs in quantum chemistry are equivalent to the problem of enumeration of perfect matchings [137]. Besides, a recent connection between dimer models and D-brane gauge theories has been discovered [54,52], providing a very powerful computational tool.
The partition function of the 2d dimer model on the square lattice was solved independently using pfaffian methods [78,41,134] for several boundary conditions, resulting in the exact calculation of correlation functions of two monomers along a row [43] or along a diagonal [55,42] in the scaling limit using Toeplitz determinants. For the general case of an arbitrary orientation, exact results are given in terms of the spin correlations of the 2d square lattice Ising model at the critical point [6,93]. Other lattice geometries have been studied as well, e.g. the triangular lattice [39], the Kagomé lattice [138,139,144], the triangular Kagomé lattice [104], the hexagonal lattice [38], the star lattice [44], or more complicated geometries [145] (see [141] for a review). The case of surface of high genius have been studied as well in [32].
The detailed study of the free energy and finite size effects began with the work of Ferdinand [40] few years after the exact solution, and has continued in a long series of articles using analytical [71,69,70,67,112] and numerical methods [95,96,94,143] for various geometries and boundary conditions. Some of these results have been motivated by the conformal interpretation of finite size effects [15,25,24] and leading to a somehow controversial result [94,124] about the central charge of the underlying field theory.
Recent advances concern the analytic solution of the problem with a single monomer on the boundary of a 2d lattice [135,142] thanks to a bijection of Temperley [133], boundary monomer correlation functions [122] and monomer localization phenomena [19,74,119]. Dimer models have regained interest because of its quantum version, the so-called quantum dimer model, originally introduced by Rokhsar and Kivelson [126] in a condensed matter context (see [108,50] for reviews) and equal to the classical model in a specific point of the parameter space. In this context, interactions between dimers have been considered at the classical and quantum level [4,3,114,33], leading to a richer phase diagram.
For the general monomer-dimer problem (cf. Appendix C for a definition) there is no exact solution except in 1d, on the complete and locally tree-like graphs [2] or scale free networks [146]. We can also mention that the matrix transfer method was used to express the partition function of the model [103,90] and a very efficient method based on variational corner transfer matrix has been found by Baxter [8], leading to precise approximations of thermodynamic quantities of the model. From a more mathematical point of view, many results exist, one can mention the location of the zeros of the partition function [58,59], series expansions of the partition function [109] and exact recursion relations [1]. Very recently an integrable version of the monomer-dimer model called monopole-dimer model has been proposed [7], sharing some qualitative features with the genuine monomer-dimer model. For d > 2 lattices, no exact solution exists for the dimer model in general, but some analytical [125,34] and numerical approaches [14,62] have been performed to study the phase diagram of the model. This lack of exact solution has been formalized in the context of computer science [75].
The prominence of the dimer model in theoretical physics and combinatorics also comes from the direct mapping between the square lattice Ising model without magnetic field and the dimer model on a decorated lattice [107,78,41,134] and conversely from the mapping of the square lattice dimer model to an eight-vertex model [9,140]. Furthermore the magnetic field Ising model can be mapped onto the general monomer-dimer model [59]. Recently, some properties of the dimer model has been proved rigorously [83][84][85] and various correlation functions have been studied as well [26,27,36]. We could also mention the study of the double dimer model [87,86] and the arctic circle phenomena [37,31] in the Aztec geometry dimer model.
Grassmann variables [12], thanks to their nilpotent properties, are very suitable to tackle combinatorial lattice models, and many of these models have been partially or entirely solved in the frame of fermionic field theory and the dimer model was one of them [127,128]. In this context, we should mention the study of spanning trees and spanning forests [22,21] as well as the edge-coloring problem [45,46] which is a special case of a more general loop model [73,92]. The approach introduced presently has been developed by Plechko in a series of papers and has been widely used to solve in a very simple and elegant fashion many problems, e.g. the 2d Ising model [117,118], boundary-field Ising model [28], the Blume-Capel model [116,30] or more general spin models [48,29].
Recently a work based on this approach has been published extending the computation of the partition function of the dimer model with an arbitrary number of monomers [5]. In this present article, we continue the analysis of this solution focusing on the effect of surfaces and corners on the free energy and correlation functions. After reviewing the classical pfaffian theory and some of its most important features in a fermionic formulation, the bosonic version of the dimer model is introduced in order to compare CFT predictions with exact calculations. A special attention will be paid to the influence of boundary conditions to the expression of the free energy and critical exponents. At the end of the article, the exact expression of monomer correlations will be employed to bring out numerical identities about the partition function of the model with boundary monomers.

Kasteleyn solution and bosonic formulation of the dimer model
In this section, the Kasteleyn-Fisher-Temperley theory of the dimer model is reminded in a fermionic field theory formulation and we show how this theory breaks down when monomers are introduced [78,134]. This theory was used to compute the partition functions, as well as dimer and monomer correlation functions in a perturbative way, leading to exact correlation exponents in the thermodynamic limit [41,43]. At the end of the section, the bosonic formulation of the dimer model is presented via the so-called height mapping. The model will be then interpreted as a Coulomb gas theory of electric and magnetic charges. This theory will be very useful in order to compare with exact results about dimer and monomer correlations performed in Section 4. Special attention will be paid to boundary conditions and corner effects in a CFT framework, which will be crucial in order to interpret finite size effects to the free energy.

Dimer model and nilpotent variables
A graph G is a pair of sets (V , E), where V is a finite set of vertices, and E is a finite set of non-oriented edges. We define the adjacent matrix (also called connectivity matrix) A = (A ij ), where the ij -entry is associated with the ordered pair of vertices (v i , v i ), then A ij = 1 if v i and v j are joined by an edge, and 0 otherwise (cf. Fig. 1 for the square lattice). The perfect matching number is the number of configurations with the property that each site of the lattice is paired with exactly one of its linked neighbors [105]. In the language of theoretical physics, the enumeration of perfect matching of a planar graph G is equivalent to computing the partition function of the dimer model on the given lattice. In the simplest form, the number of dimers is the same in all the configurations, and the partition function is given by the equally-weighted average over all possible dimer configurations. 1 In the following, we will include equal fugacities t for dimers, so that the average to be taken then includes weighting factors for dimers and we write the partition function as where the Hamiltonian for the dimer written using commuting nilpotent variables (see Appendix A) can be written as a sum over every vertices (see Fig. 2), preventing two dimers to occupy the same site where A ij is the adjacent matrix of the lattice considered. Let us put β = 1 in the following. The nilpotent variables can be seen as commuting Grassmann variables, or simply a product of two sets of standard Grassmann variables where η i = θ iθi . The perfect matching number of the graph G is equal to the partition function in the case t = 1 Fig. 2. At every vertices, we put a nilpotent variable η, such that η 2 = 0 forbidden two dimers to occupy the same site.
We report the reader to Appendix A for the definition of the haffnian of a matrix. In the second line we decomposed the nilpotent variables using two sets of Grassmann variables, and we finally found the well known graph theory result perfect G = hfA.
Considering holes in the perfect matching problem is equivalent to removing rows and columns at the positions of the holes in the adjacent matrix (see Fig. 3). The resulting combinatorial problem is called the near-perfect matching problem. The partition function of the dimer model with a fixed number of holes (monomers) can be written as where the index n stands for the number of monomers in Q n [t]. Finally, the result is the same but now the matrix A \{q p } is the adjacent matrix of the original graph with positions of the n monomers removed. Suppose we remove two sites q 1 and q 2 on the graph G, then it is similar to introduce two nilpotent variables η q 1 and η q 2 on the lattice, the correlation function between these two monomers is then and more generally the n-point correlation function reads The partition function and correlations can be studied in the case t = 1 as well, in that case, the matrix elements of A are a ij = ±t and the generalization is straightforward. Generally correlations between monomers are equal to correlations between nilpotent variables in this framework, which can be written in terms of a ratio between two haffnian. Unlike the determinant which can be computed by an O(L 3 ) time algorithm by Gauss elimination, there is no polynomial time algorithm for computing permanent. The problem of converting a permanent problem into a determinant problem is a long standing problem in pure mathematics, the simplest version of this problem is called the Pólya permanent problem [120]. Given a (0, 1)-matrix A := (A ij ) L×L , can we find a matrix B := (B ij ) L×L such that permA = detB (or equivalently hfA = pfB) where B ij = ±A ij .

Haffnian to Pfaffian conversion and Kasteleyn solution
The close-packed dimer model can be solved on any planar lattice by using Pfaffian techniques. These techniques were introduced in the early sixties by Kasteleyn [78] and Temperley [134] and give an answer to the Pólya permanent problem for some particular conditions. The Kasteleyn theorem is a recipe to find a matrix 2 K in such way that hfA = pfK, where the elements of the matrix are K ij = ±1. Kasteleyn theorem is based on a special disposition of arrows on the edges of a planar graph, 3 the product of their orientation around any even-length closed path should be −1. Such a disposition is given in Fig. 4(b) for the square lattice. We define an antisymmetric matrix K, where if the arrow points from i to j −1 if the arrow points from j to i 0 otherwise Kasteleyn theorem states that the perfect matching number of a given planar graph G is given by which is equal to ± √ det K. The ± sign is chosen to make the perfect matching number positive, henceforth we will omit this sign in the rest of the article. Differently, this pfaffian can be express in terms of Grassmann variables (cf. Appendix A) For the square lattice, we can choose Boltzmann weights t x and t y for horizontal and vertical dimers, then the pfaffian can be computed by using Fourier transform and Kasteleyn found for free boundary conditions (cf. [78] for details on calculations) In Table 1, we compute Q 0 [1, 1] using MATHEMATICA ® , for different M and N with t x = t y = 1 (perfect matching number). All these values can be numerically checked using diverse algorithms which enumerate all the possible configurations on the square lattice (see [100] for details). In the rest of the paper we will omit the labels t x and t y in the partition function and just keep Q 0 .

Entropy in the thermodynamic limit
The asymptotic form L → ∞ of the partition function (for M = N = L) can be easily found from Eq. (10) where G is the Catalan constant. 4 The factor G/π is the entropy by site of the dimer model on the square lattice. This entropy can also be calculated for other bipartite lattice like the honeycomb lattice, and for non-bipartite lattice like triangular, Kagomé lattice and triangular Kagomé lattice. The pfaffian method can be used to compute the partition function of the dimer model on various geometries and boundary conditions leading to different values of the entropy (cf. [141] for review), as long as the lattice has a Kasteleyn orientation (i.e. planar according to the Kasteleyn theorem) and as long as N × M is even. Obviously it is impossible to fill an odd size lattice with dimers, without leaving one site empty. We shall take notice later that the form of the free energy is strongly dependent on the parity of the lattice.

Probabilities and correlations
The Kasteleyn matrix is a very powerful object, which gives us all the details about probabilities of presence of dimers [43]. For example the occupation probability P[i → j ] of a dimer on the link ij is where K ij and K −1 ij are the corresponding matrix elements of the Kasteleyn matrix and its inverse. The probability of two dimers on the links ij and mn is In the rest of the paper, we shall use the term correlation even though this quantity is a normalized probability. Correlations between more than two dimers are available using the Kasteleyn matrix as well. It can be shown [41,43] that dimer-dimer correlations on the square lattice decrease as the inverse square of the distance between the two dimers in the thermodynamic limit Furthermore it has been also shown that correlations are always critical [79] for bipartite lattices and a contrario exponential for non-bipartite lattices as the triangular lattice for example, where the fermionic theory underlying is a massive theory [39].

Monomer correlation functions
Throughout this work the monomer-monomer correlation function C will be defined as the ratio of the number of configurations with monomers at fixed positions to the number of configurations without monomers. Thus computing a monomer-monomer correlation is stricto sensu equivalent to compute the partition function with two sites (and all the links connected to these sites) deleted. Since such a graph is still planar, Kasteleyn's construction is still applicable. The one complication is that we must ensure that on the new lattice with deleted sites, the number of arrows is still clockwise-odd. If all the monomers are located on the boundary of the lattice at ordinate {x i }, there is no non-local defect lines between monomers (see Fig. 5), and the modified matrix K \{x i } defined from K by removing all the rows and columns corresponding to the monomer positions has still the proper Kasteleyn orientation. Then the pfaffian of this modified Kasteleyn matrix K \{x i } gives us the partition function of the dimer model with fixed monomer positions. It follows that the correlation function between two monomers on the boundary is where x 1 and x 2 are the positions of the two monomers. This pfaffian has been computed by Priezzhev and Ruelle (cf. [122] for details) in the thermodynamic limit for an arbitrary number of monomers at positions {x i }, using a perturbative analysis of the matrix K \{x i } around the original Kasteleyn matrix K. The result for the 2n-point correlation is given by where the matrix element C ij := C(x i , x j ) is the 2-point function of a 1d complex free-fermion, equal to if x i and x j are on opposite sublattices and C ij = 0 otherwise. For monomers in the bulk, the things are much more complicated. One sees that the product of arrows around a deleted site is now equal to +1 (see Fig. 5). We thus must construct a string of reversed arrows from one monomer to the second (see Fig. 5(b)). As long as the arrows are chosen to make all plaquette clockwise odd, the correlation is independent of the choice of the path. In the general case of bulk monomers, the relation (15) is no longer correct, because the matrix K \(x i ,x j ) is no more a Kasteleyn matrix. Then correlations betweens two monomers defined by Q 2 (x i , x j )/Q 0 are not equal to correlations between two Grassmann variables a i a j , but disorder operators must be added where the sum is over all the links connecting sites i and j , to take account of the reversing line between the two monomers. Using a pfaffian perturbative analysis it was shown that monomer correlations decrease at the thermodynamic limit as [41,43] This result is very similar to the construction of the spin correlation functions in the Ising model in terms of fermionic variables [77,121]. In fact, on the square lattice, the monomer-monomer correlations were shown to have the same long-distance behavior as the spin-spin correlations in two decoupled Ising models, explained by the deep relation between the two correlation functions for the square lattice given by Perk and Au-Yang [6]. These disorder operators are absent in the haffnian theory, Eq. (7), and are the price to pay to solve the problem analytically.

Height mapping and Coulomb gas formalism
To any dimer covering we can associate a height on the dual lattice (on the plaquette) which is defined as follows [147,102,64,61]. When encircling an even vertex in the positive (counterclockwise) direction, the height h increases by +1 upon crossing an empty edge and decreases by −3 upon crossing an edge that is covered by a dimer (cf. Fig. 6). It is easy to notice that for the allowed configuration the average values h vertex that the height variables can take at a given site of the direct lattice (a vertex) are h vertex = ±3/2 or h vertex = ±1/2. 5 On the other hand, a uniform shift of all the heights by one unit leads to an equivalent state. This mapping works stricto sensu for the close packed case. We will find it simpler to work with the rescaled height field φ = π 2 h. By fixing the rescaled height at an arbitrary point, e.g. φ(0) = 0, these rules define the entire height function φ( r) uniquely. By integrating out the short distance fluctuations, one obtains an effective quadratic action S for the bulk height field φ( r), defined in the continuum, which corresponds to the long-wavelength modes Here g a constant which controls the stiffness of the height model. It is a priori unknown. The field has to be invariant under the transformation φ = φ + 2πn to respect lattice symmetries.
The derivation of this Gaussian field has been actually done rigorously by mathematicians [84,85]. Electric charges e correspond to vertex operators appearing in the Fourier expansion of any operator periodic in the height field. Dual magnetic charges m correspond to a dislocation in the height field and correspond to the dual vertex operator V e (z) = : e ieφ : where ψ is the dual field of φ and defined as These operators are primary operators of the c = 1 CFT. The scaling dimension associated to the insertion of a particle with electromagnetic charge (e, m) is given by For example, two monomers on opposite sublattices correspond to two charges m = 1 and m = −1. It is known from exact results [41,43] that, the exponent for bulk monomer-monomer correlations is 1/2, it fixes g free = 1/4π for the stiffness constant of the Gaussian field theory describing the free dimer model. We saw previously that bulk dimer-dimer exponent is 2. Hence the bulk monomer and dimer scaling dimensions defined by x In this theory the conformal spin of an operator is defined by s(e, m) = em, then monomers and dimers are spinless particles but fermions which are order-disorder composite operators [77] have magnetic and electric charges and carry spins 1/2. The fermion operator has then scaling dimension (1/2, 1) = 1/2. It is also possible to define parafermion operators which obey fractional statistics [51] for particular values of the stiffness g. The use of such a mapping to study correlation functions dates back to Blöte, Hilhorst, and Nienhuis [111,16]. The neutral 2-point correlation functions for vertex operators are then given by the standard formula [110] The general monomer 2n-point function is given by the product where {b k } and {w k } are the sets of the even/odd sublattice coordinates. In this Coulomb gas interpretation of the dimer model, monomers located on the same sublattice are seen as repealing equal charges, and monomers on opposite sublattice are seen as attractive opposite charges. It is known from exact results [122,5] that, the exponent for monomer correlations is 1 in the case (see Eq. (16)) where the monomers are restricted to live on a boundary. Hence the surface monomer scaling dimension is x (m) s := 1/2.

Rectangular geometry and boundary conditions
On Fig. 7, we show a configuration of dimers on the square lattice with free boundary conditions. The choice of these conditions imposes very specific boundary conditions for the height variable. Indeed, we observe that the value of the height h is (..101..) along the boundary until the corner and (.. − 10 − 1..) along the next boundary [129] (a similar analysis has been done for the strip geometry [131,132]). Then we have, from the point of view of the field theory two different averaged boundary conditions either sides of the corner, let us call the boundary fields h b = 1/2 and h b = −1/2. The proper way to understand how corners change the behavior of the height field is through boundary CFT (cf. [60] for introduction). Hence it is natural to introduce local operators [23] acting on the corners of the domain, these abjects called boundary condition changing operators (bcc) can be showed to be primary operators of the CFT [24]. A careful analysis shows that the difference of boundary conditions has non-negligible consequences in the thermodynamic limit and the precise procedure based on the contribution of theses operators on the free energy can be found in [131]. It has been shown that the dimension of the bcc operator which creates a field shift of value with g free = 1/4π here. In the previous situation without any monomer, the corner shift in the height h is equal to h b = 1 then φ b = π/2, hence the dimension of the corner bcc operator is h bcc = 1/32. Furthermore, the addition of a monomer on the boundary or at the corner will change the value of the field, and we will behold further that it will be relevant in order to study quantities as the free energy and correlation functions. A general framework to study partition functions and conformal boundary states on the rectangular geometry with different boundary conditions in a boundary CFT framework has been developed recently [17,18]. In Section 4, finite size effects to the free energy for the free and interacting cases will be study in this CFT framework, and the influence of these bcc operators will be crucial to identify the correct underlying central charge of the theory. Furthermore the presence of monomers on the boundaries or at the corners will change again the values of the bcc operators which will be important for the study of surface and corner correlation functions.

Exact partition function and corner free energy
In this section, the fully-detailed Grassmann solution of the dimer model with an arbitrary number of monomers is presented, which will lead to the exact form of correlation functions in a pfaffian formulation. In particular, we show that the problem become simpler when we consider boundary monomers, and a closed expression for correlations can be found. Hereinafter the solution of the dimer model on an odd size lattice with one boundary monomer is introduced, in agreement with the Tzeng-Wu solution [135]. Those solutions will be used to extract finite-size scaling behaviors of the free energy in a CFT framework, where boundary changing conditions operators have to be carefully studied to infer the central charge of the model, in contradiction with a recent article [68].

Plechko pfaffian solution
Here we just recall the framework of the Plechko solution [56,57] of the dimer model. As we have seen in Section 2, the partition function can be written for a general graph using nilpotent variables We are now working on the square lattice with free boundary conditions, then the partition function reads where t x and t y are the horizontal and vertical Boltzmann weight, and m and n refer to the coordinates. The partition function can be written using Grassmann variables (see Appendix B for details), this leads to a block representation of the action in the momentum space, for momenta inside the reduced sector 1 ≤ p, q ≤ L/2. The four components of these vectors will be written c μ α with μ = 1 · · · 4, leading to 6 where the antisymmetric matrix M is defined by This matrix can be written as where the matrices x and y are (119) is directly related to the pfaffian of the matrix M (cf. Appendix A for details) and is simply equal to Finally one simply obtains the following well known result The fermionization can also be performed for toroidal boundary conditions. We refer here to the experience with the 2d Ising model on a torus [117]. The final result can be written in terms of a combination of the periodic-antiperiodic boundary conditions for fermions c M+1n=±c 1n and c mN+1=±c m1 .

General case
Let us now consider the case where an even 7 number of monomers are present on the lattice at different positions r i = (m i , n i ) with i = 1, · · · , 2n. The partition function Q 2n ({r i }) is the number of all possible configurations with the constraint imposed by the fixed monomers. This quantity can be evaluated by inserting nilpotent variables η m i n i in the partition function, which prevents possible dimers to occupy sites r i . It can be useful to introduce an additional Grassmann variable h i such that η m i n i = dh i exp(h i η m i n i ). These insertions are performed at point r i in Q 0 , and the integration over η m i n i modifies L m i n i → L m i n i + h i . However, by moving the dh i variable to the left of the remaining ordered product, a minus sign is introduced in front of each b mn i −1 in B mn i for all m > m i . We can replace b mn−1 by mnbmn−1 such that mn i = −1 for m > m i , and mn = 1 otherwise. The integration is then performed on the remaining variables (a, ā, b, b ) as usual, so that Q 2n ({r i }) can be expressed as a Gaussian form, with a sum of terms corresponding to the monomer insertion, resulting to the following form for the partition function where D[c]D[c] = mn dc mn i dh i where the index i runs over the set of positions r i = (m i , n i ) of the 2n monomers. The inclusion or monomers is equivalent to inserting a magnetic field h i at points r i , as well as a sum of quadratic terms c mn i −1 c mn i running from the hole position to the boundary on the right (see Fig. 8). Another possibility would be to join two monomers by a line of terms by moving dh m i n i until dh m j n j as in the Kasteleyn theory. In this case, the additional quadratic terms in the action starting from r i and ending on the boundary have to be treated in the computation of the Grassmannian integral. We first rewrite S 0 in the Fourier space using the block partition label α = (p, q) for momenta p and q inside the reduced sector 1 · · · L/2, and vectors c α = t (c pq , c −pq , c p−q , c −p−q ). Also the 4 components of vector c α will be written c where the antisymmetric quadratic form M α is defined by Eq. (31). The part of the field interaction can be Fourier transform as before with a linear field H pq depending on h i s and we obtain The last contribution connecting the monomers to the boundary can be written as i 2 c μ α V μν αβ c ν β , with matrix V αβ = V pq,p q given by The different components V μν αβ are given implicitly, for the first elements, by V 11 αβ = V pq,p q , V 12 αβ = V pq,−p q , V 21 αβ = V −pq,p q , and so on. Then the total fermionic action contains three terms The first two terms contain only modes of the same sector α, and the last connects modes from different sectors α and β. Matrices  Fig. 9(a)). By construction, W is antisymmetric and satisfies W μν αβ = −W νμ βα . This matrix can be represented as a block matrix of global size (1,2) V (1,1), (1,3) · · · −V (1,1), (1,2) M (1,2) V (1,2), (1,3) · · · −V (1,1), (1,3) −V (1,2), (1,3) M (1,3) · · · . . .
where each of the (L 2 /4) × (L 2 /4) blocks is a 4 × 4 matrix. Labels α are ordered with increasing momentum (1, 1), (1, 2) · · · (1, L/2), (2, 1) · · ·. In the action the linear terms in c μ α can be removed using a linear change of variables c μ α → c μ α + g μ α , where g μ α are Grassmann constants. After some algebra, we find the correct values for the constants that eliminate the linear contribution are given by g After substitution of these values in the overall integral, and a rescaling of variables c α → c α / √ i, the action becomes a product over variables c and h i can be rewritten using a 4-dimensional vector, such as where the vector depends only on the location parity of the monomer r i in the bulk. Functions r i (α) are normalized α r i (α)r j (α) = δ ij . Then we obtain the following formal and compact expression for Q 2n ({r i }) Finally, we found a pfaffian expression of the partition function We verify easily that the matrix C is antisymmetric by using the antisymmetry property of W or Then, C ij can be formally expressed as a scalar product C ij = α,β i,α |W −1 α,β | j,β . Q 2n is therefore a product of two pfaffians where the monomer locations are specified in matrix W . The matrix V can be rewritten using additional matrices after considering the different components (μ, ν). We can indeed express V using four functions u a,s k (α, β), and v a,s k (α, β), for each monomer at location r k = (m k , n k ), with m k < L, and such that where the cc are defined as Functions u and v are given by More explicitly the result is This close solution to the dimer model with an arbitrary number of monomers at fixed location can be formally used to get access to some informations about the general monomer-dimer model (see Appendix C).

Boundary monomers
Although the general problem remains in principle tractable, we can simplify it further by considering monomers on the boundary m i = L of the rectangle only, with the convention n i > n j if i > j. When the monomers are on the boundaries of the domain, the off diagonal term V vanishes and the matrix W reduces to the matrix M (cf. Fig. 9(b)). In that case, the previous action contains no defect line term, only the couplings with fields remain and the problem can be transformed into a 1d system of particles on a chain (see Fig. 10). Using the previous Fourier transform for the c's variables, the magnetic field term becomes Then the partition function can be integrated on the c mn 's using the block partition (p, q) and we find that where the prime symbol is meant for summation over half of the modes p, q = 1 · · · L/2. This action would vanish if all the n i were for example even, since the Grassmann fields satisfy in this case H pq H p−q = −H 2 pq = 0. In general, the field action can be rewritten as a quadratic form over the real-space fields h i : S H = i<j C ij h i h j , where the elements of the matrix C are antisymmetric C ij = −C ji , and equal to Fig. 11. Diagrammatical representation of the pfaffian of C for 2, 4 and 6 monomers.
These elements are zero if n i and n j have the same parity and in general the integration over the field variables h i leads directly to a pfaffian form for the partition function Q 2n ({r i }) = Q 0 (−1) n pf(C). The (−1) n factor comes from the rearrangement of the measure dh 1 · · · dh 2n = (−1) n dh 2n · · · dh 1 , so that the 2n-function reads This sign could also be absorbed in the definition of matrix elements C ij → −C ij . The resulting partition function is always positive with this definition. For example, if there are 2 monomers at the boundary, pf(Ĉ) = −C 12 , and 4 monomers 8 leads to pf(Ĉ) = C 12 C 34 − C 13 C 24 + C 14 C 23 (cf. Fig. 11). We should notice that the expression (53) is an exact closed expression for the 2-point correlation between monomers on the boundary m i = m j = L, leading to an explicit result for the partition function with 2, 4...2n monomers. We will show in a next section how to find the correlation between monomers on different (opposite or adjacent) boundaries of the lattice.

Single monomer on the boundary
We can recover the partition function for one monomer on the boundary using the previous analysis. In this case the size L has to be odd in order to accommodate for the presence of one single monomer. The action (119) is still valid, but the Fourier transform leads to a different block arrangement for the bulk terms in Eq. (123) which are represented by the red zones in Fig. 12 S 0 = 2it x  S 0 contains Fourier modes that cover the Brillouin zone (see Fig. 12) except for term c 1 2 (L+1) 1 2 (L+1) which is located in the middle of the zone and not present in the sums (55). The integration over Grassmann variables c pq is therefore zero in absence of coupling with this mode. Inserting a single monomer on the boundary at location r = (L, n) is equivalent, as previously demonstrated, by inserting a Grassmann field h and a field contribution S H = c Ln h in the action, which can be expanded using Fourier transformation Since there is only one Grassmann field, all terms H pq are only proportional to h, and therefore the quadratic form in Eq. (52)  2 cos( πp L+1 ) = L+1 2 , to simplify the previous expression and recover the Tzeng-Wu [135] solution We can notice that the partition function with one monomer in a system of size L × L (L odd) is equal to the partition function without monomers on a lattice of size L − 1 × L − 1. The probability is therefore constant for all location of the monomer, at even sites only, the last term in bracket being equal to zero (n odd) or unity (n even), proving that the monomer is fully delocalized on the boundary, unlike the bulk case where monomers are actually localized in a finite region of the domain [19,74,119].

Corner free energy and the central charge controversy
The study of finite size effects in statistical physics is a long standing and still active field of research [123]. A fortiori the possibility to solve a model in a non-homogeneous geometry [63] is of prime interest for the understanding of behavior of physical systems in real situations. In the case of the dimer model on the rectangle with free boundary conditions, the system admits surfaces and corners, both of them play an important role in the behavior of the free energy in the thermodynamic limit. The exact solution of the close packing dimer model equation (10) on an even lattice (M × N = 2p) allows for the study of the finite size effect of the free energy, and the finite size analysis has already been performed in the early time of the dimer model history [40].
Furthermore the exact solution (56) of the dimer model on an odd lattice (M × N = 2p + 1) with a monomer on a boundary allows for the study of the free energy in that case as well. Adding a finite number of monomers in the dimer model is equivalent to a zero density of monomers in the continuum limit. Hence the presence of monomers does not give any contribution in the expression of the free energy, the only feature which plays a crucial role is the parity of the size of the lattice (even or odd). Because they are the simplest expressions of an even (odd) lattice with an even (odd) number of monomers, these two partition functions are sufficient to study all the details of the asymptotic limit of the free energy. In the following let us choose the square geometry M = N = L for simplicity. 9 The free energy on a finite lattice of typical length L at criticality has the generic following form 10 The first term is the extensive contribution of the free energy, whereas the second one represents contribution from the lattice surface. In general, the coefficients f bulk and f surface are non-universal, but the coefficient f 0 is assumed to be universal, depending only on the shape and boundary conditions of the system. Universal properties of critical models appear in the subleading corrections and the value of f 0 is known to be simply related to the central charge c of the underlying conformal field theory. The study of statistical systems and their field theory representation in the presence of corners has been covered extensively, e.g. Ising and Potts model, loop model and percolation [72,136,99,98], using various theoretical and numerical machineries. In two dimensions, as pointed out by Cardy and Peschel [25], the universal contribution to the free energy of a critical system in a domain with a corner with angle θ has been determined using the complex transformation which maps the upper half-plane onto the corner and looking at the holomorphic component of the stress-energy tensor in the corner. This mapping gives us the explicit form of the logarithmic contribution F corner = a log L in Eq. (57) and the result is c 24 θ π − π θ log L. It turns out that an additional complication arises because of the bcc operators [24] acting in the corners (cf. Fig. 13) as we saw in Section 2. In that particular case, the Cardy-Peschel contribution is slightly modified to taking into account this change of boundary conditions [89], and the logarithmic corner contribution becomes 9 The entire procedure can be extended to the case M = N where the aspect ratio has to be taken into account, and no significant change appears [88]. 10 The presentation here closely follows [130]. where h bcc is the scaling dimension of the bcc operator. This bcc operator changes the boundary condition either sides of the corner in the height field representation. In their recent paper about corner free energy contribution in the free boundary conditions dimer model with one monomer at the boundary [68], the authors analyzed the asymptotic contribution of the four corners of the rectangular system without taking into account this bcc operator, i.e. taking h bcc = 0 in the previous formula, and concluded that the central charge of the dimer model is c = −2. It is beyond the scope of this paper to enter into details in the wide area of conformal analysis of finite size effects (see [130] for more details), nor all the literature of c = −2 models, but just to pointing out the difference of result when one looks carefully at the bcc operator contribution in the corner free energy. In their paper, the authors found that the contribution of the four corners of the L × L (L odd) lattice is 11 The CFT formula (59) gives in the square geometry case (θ = π/2) where h (ν=1..4) bcc are the dimensions of the four bcc operators living on the corners. Taking h bcc = 0 leads de facto to c = −2, suggesting that the dimer model may be a logarithmic CFT (LCFT) [71,122,124]. This statement is also based on the mapping of the dimer model to the spanning tree model [133] and, equivalently, to the Abelian sandpile model [106] which both belong to a c = −2 LCFT, facts which we do not dispute here. The problem with this analysis is the oversight of the bcc operators acting on the corners. We know that the partition function with one monomer on the boundary does not depend on the location of the monomer, then let us choose to put it on the corner for simplicity. In this following case, the height field is shifted either sides of the corner, precisely we obtain the value h = (..0101..) in one side and h = (.. − 2 − 3 − 2 − 3..) on the other side (see Fig. 14), making a height shift of φ b = 3π/2, therefore using Eq. (27), 11 Cf. [68] for details of the asymptotic calculation. the dimension of the bcc operator on the corner is 9/32 (which is actually the dimension of the original bcc plus the dimension of the corner monomer operator as we will see in the next section). The three other corners induce a height shift of π/2, thus the dimension of the bcc operators is 1/32 as in the case where monomer are absent. Finally we find The comparison with the asymptotic result (60) gives us the value of the central charge of the free bosonic field theory, i.e. c = 1. 12 Obviously this result seems completely natural in the height mapping framework or in the free complex fermion representation of the dimer model equation (37), nonetheless the presence of bcc operators acting on the corners has never been extensively studied in this context, leading to misinterpretation of asymptotic results and thus to a different value of the central charge of the theory. In the pure dimer situation Eq. (10), the four bcc operators has dimension 1/32, we claim that the corner free energy should be equal to This result nicely agrees with the literature, where no corner free energy term has never been found in the free boundary conditions close packed dimer model, neither theoretically [20] nor numerically [94], strengthening our analysis. 13 Previously another type of finite size effect analysis [71] has been performed for the close packing dimer model on a strip with periodic and free boundary conditions, and it has been shown that the result depends strongly on the parity of the length of the strip as it should be, and the notion of effective central charge has to be introduce [66]. Though the comparison with the CFT result gives the value c = −2, we claim here that the analysis of the free energy in the height mapping formulation of boundary conditions, leads again to c = 1. In the following section, we shall compute analytically correlation functions between monomers where boundary fields has to be properly interpreted and we will show that the exact solution fully agrees with this theory, reinforcing the present result about the significance of these bcc operators. Finally we conclude, that the free bosonic formulation of the dimer model allows for the complete study and interpretation of finite size effects in a CFT context, and unified all the results known about this very simple but not trivial model.

Exact correlations: discrete and continuous cases
In this section, detailed computations of correlations are performed in terms of disorder operators. A particular attention will be paid to the special case of boundary monomers, where a closed-form expression is obtained, valid for any of the four boundaries of the rectangle. Then, numerical evaluations of the exact pfaffian solution are done and diverse monomer and dimer 12 One can notice that the same result holds if the monomer is somewhere on the boundary but not at the corner, in that case, there is five bcc operators, four on each corners of dimension 1/32 and another responsible for the shift at the surface with dimension 1/2 then using Eq. (59) in accordance with the fact that the partition function with one boundary monomer does not depend on its location. 13 The same analysis can be done for the same model with periodic boundary conditions in one direction and free in the other, leading to the same conclusion about the central charge.
correlations are obtained for bulk, surface and corner cases, leading to the whole set of scaling dimensions which can be compared to the bosonic theory with g = 1/4π .

Fermion correlations and disorder operators
The addition of monomers in the dimer model is therefore equivalent to inserting a magnetic field h i at points r i , as well as a line of defect running from the monomer position to the right boundary m = L. If two monomers have the same ordinate n i = n j , the line of defects will only run between the two monomers and will not reach the boundary. This can be viewed as an operator acting on the links crossed by the line and running from a point on the dual lattice to the boundary on the right-hand side. More specifically, we can express the correlation functions, after integration over the fermionic magnetic fields h i , as an average over composite fields where μ(r + e 4 ) is a fermionic disorder operator, whose role is to change the sign of the vertical links across its path starting from vector r + e 4 on the dual lattice (cf. [5] for details). The integration · · · 0 is performed relatively to the action S 0 . Likewise the Kasteleyn theory, where disorder lines are absent on the boundary and where correlations between monomers correspond to correlations between Grassmann variables, Eq. (15), here the correlation between monomers on the boundaries are exactly correlation functions between the fermionic fields This result about monomer correlations written in terms of disorder operators is the fermionized version of the Coulomb gas framework, where monomers act like dual magnetic charges which create a dislocation of the height field and correspond to the vertex operator of the corresponding bosonic field theory.

Perturbative expansion of the 2-point function
In the case where the monomers are on the boundaries, we were able to compute exactly the 2-point correlation function (53) in the discrete case. In the bulk case, the things are much more complicated and an exact closed-form expression on the discrete level seems out of reach. Nevertheless, a perturbative expansion can be performed to evaluate the pfaffian expression of the correlation function. We start from the exact pfaffian expression of matrix C The inverse matrix W −1 can be computed using formally the expansion In particular it is convenient to write the inverse matrix M −1 as withā In the following we will consider only the first of this expansion The structure of this expansion make possible a further diagrammatic expansion of the quantities C ij as a series of term C (2) ij := c ij γ (2) ij . (72) It is easy to see that i | x | j = 0 for monomers on the boundary or on the same column, when m i = m j . A contrario, for pairs of monomers on the same line n i = n j , we have i | y | j = 0.
The first term C (0) ij can be expressed in the discrete case as This expression is valid for 2 monomers on any of the four boundaries of the lattice, and is identical, when m i = m j = L, to expression obtained for the same (m i = n i = L) boundary case, Eq. (53). Indeed the first order of the expansion (68) is valid only on the boundaries where the matrix W is actually equal to the matrix M. One could demonstrate that this 2-point correlation is actually exact C (0) ij = C ij for boundary monomers because of the cancelation of higher terms in the perturbative expansion, accordingly the expression (73) is a general exact result for any positions anywhere on the four boundaries. Therefore, it will be very efficient to use this exact closed-form to evaluate scaling behaviors of correlation functions between monomers on the surface and at the corners. This present perturbative expansion can be performed to the next leading order to evaluate bulk correlations, and will be detailed elsewhere.

Scaling behavior of monomer correlation functions
Here, we shall analyze monomer-monomer correlation functions using our pfaffian solution detailed previously and compared to the Coulomb gas interpretation of the dimer model. As we saw in Section 2, the dimer model on a rectangular geometry admits a bcc operator on every of the four corners, and it has to be taken into account for the analysis of the scaling dimension operators, in particular for corner correlations. Indeed in the case of monomers deep in the bulk 14 or deep in the surface, 15 Fig. 16(a)) These known results are in perfect agreement with our exact solution (cf. Fig. 15) where we fixed the positions of two monomers for increasing system size L (all the correlations are measured for t x = t y = 1). The behaviors of bulk and surface monomer correlation functions had already been studied in several papers, and the scaling dimensions are related to the scaling dimensions of operators of the 2d Ising model via the expression of monomer-monomer correlations as spinspin correlations [6]. At present, we consider the effects of corners in our system, which seems to be more difficult to consider as we have seen for the corner contribution to the free energy. Fortunately conformal invariance predicts a relation for the scaling dimension of an operator in the vicinity of a corner of an angle θ in terms of the scaling dimension of the same operator on the surface [25] x c = π θ x s .
If we believe in this formula, we should obtain the value x (m) c = 1 for θ = π/2, which leads to the behavior C(L) ∼ L −2 for corner-corner correlation functions. This result contradicts our exact evaluation, where the exponent seems to change according to the exact location of the monomers (see Fig. 16(b) and Fig. 18), and where three different cases arise. Unlike the surface and bulk cases where the scaling dimensions are uniquely defined, the corner scaling dimension appears to be less trivial to analyze, and the influence of the bcc operators has to be carefully taking into account. We should mention that the same kind of analysis has been done for the Ising model, where the magnetization was measured for various spins close to a corner [115]. We saw previously in Eq. (61), that bcc operators add a logarithm term in the expression of the free energy F of a rectangular system, then this contribution to the partition function scale as Indeed putting two monomers exactly on the corner of the same boundary (see Fig. 18(a)) is equivalent to a height shift of value 3π/2 in each corner. This height shift is induced by an operator of dimension 9/32, leading to the following behavior of the partition function with the two monomers Here the correlation function scale then as C(L) = Q 2 Q −1 0 ∼ L −1 leading to the value x (m) c = 1/2 of the monomer corner scaling dimension. Nevertheless, if one choose the diagonal corners (see Fig. 18(b)), we place a monomer on the first corner which is again equivalent to the insertion of a bcc operator of dimension 9/32 (see Fig. 14) and the other one on a neighboring site of the other corner which is equivalent to the insertion of a bcc operator of dimension 25/32 (see Fig. 17), we found the behavior C(L) ∼ L −2 . Finally, if both of the two monomers are on a neighboring site of a corner (see Fig. 18(c)), then the correlation function is C(L) ∼ L −3 . This three different situations are summarized in Fig. 18, showing that our exact computations are  in perfect agreements with Coulomb gas predictions for the behavior of correlation functions in the vicinity of a corner. A more general statement is that the monomer scaling dimension near a corner depends crucially on the sublattice considered as explained in Fig. 16(b). This phenomena leads to two different values of the scaling dimension for corner monomers x (m) c = 1/2 or 3/2 which is in agreement with the CFT formula (75) in average when lattice effects are forgotten. A general study of the finite size behavior of correlation functions can be performed as well, leading to the following scaling ansatz where r 1 and r 2 are the positions of the two monomers, with respective scaling dimensions x 1 and x 2 . The scaling function (u) depends on the position of the operators and goes to a constant in  the scaling limit u → ∞ (see Fig. 19). The translation and rotational invariance has been checked analytically and numerically, and in the following we will use |r 1 − r 2 | = r This scaling behavior is shown for bulk-bulk and surface-surface correlations in Fig. 19. The exact form of the scaling function seems hard to obtain explicitly, but at least for boundary and corner cases it should be possible to extract the scaling behaviors using the expression (73) of exact correlation functions. We let this question for a future work and we hope that comparisons with CFT predictions can be made.

Scaling behavior of dimer correlation functions
As we have shown in Section 2 of this article, the bulk correlation between a dimer covering the two neighboring sites i and j and another dimer covering the two neighboring sites m and n can be computed in the Kasteleyn-Fisher-Temperley pfaffian formalism leading to a behavior in L −2 in the thermodynamic limit. Let us note D(L) this quantity. In the Coulomb gas approach, dimers are interpreted as electric charges with scaling dimensions x Actually it is straightforward to study dimer correlations here. Indeed a dimer can be seen as two neighboring monomers, thus a dimer-dimer correlation is simply a 4-point monomer-monomer correlation, which can be evaluated with our solution. In the following, one shows how to construct the dimer-dimer correlation in the boundary case, for the bulk case the situation is essentially the same but expressions are less convenient. Explicitly the correlation between two dimers at position (r i , r j ) and (r m , r n ) is If we choose that r i and r m are on the same sublattice (then r j and r n are on the other one), a straightforward consequence is that the second term C im C jn vanishes, moreover the first term C ij C mn tends to a constant in the thermodynamic limit in such a way that we can define the dimer-dimer correlation function as In this way, all the configurations of dimer correlations are available, and all the scaling dimensions of electric charges (bulk, surface, corner) may be analyzed stricto sensu and compared with the Coulomb gas theory. We can show that the well known bulk behavior is recovered very precisely, furthermore, surface and corner correlations may be examined as well leading to the following behaviors Unlike monomer correlations, dimer correlations are much easier to interpret in the Coulomb gas framework. Indeed, the absence of additional change of boundary conditions 16 in the partition function allows for a direct determination of dimer scaling dimensions. The particular form of dimer correlations (81) predicts that x (d) 18 We notice here that the formula (75) checked out in that case. A careful and detailed study of surface and corner operators has to be performed to unravel this point. The scaling form of correlation functions (78) holds in the dimer case as well, using the dimer scaling dimensions (see Table 2). The solution presented in this article can be also used to calculate more complex correlation 16 Of course the four corner bcc's with dimension 1/32 are still present but do not play any role in dimer-dimer correlation functions. 17 Let us notice that the fact that x s is a pure coincidence. 18 Table 2 Bulk, surface and corner values of dimer and monomer scaling dimensions for the free (g free = 1/4π ) fixed point. The corner monomer scaling dimension depends on its exact location.
Scaling dimension (g free = 1/4π ) Bulk Surface Corner functions, combining dimer and monomer scaling dimensions. A posteriori, more complicated object like trimers, quadrimers or more generally, string of k neighboring monomers (k-mer) can be studied as well, which correspond to various charged particles in the Coulomb gas formalism.

About some combinatorial properties
In this section, one shows a curious combinatorial analogy between the partition function of the close packing dimer model on an L × L square lattice with open boundary conditions, and the same partition function with boundary monomers. One start reminding some properties of the pure dimer model partition function, and we show, thanks to our exact calculation of the partition function with 2n monomers, that this analogy can be understood and demonstrated. Hereinafter, the Boltzmann weights t x and t y are taken to be the unity in such way that the partition function is exactly equal to the perfect matching number. All the results presented in this section have been checked with depth-first [100] algorithms up to size L = 10. For bigger sizes, Monte-Carlo simulations [101] or transfer matrix calculation [103] has to be implemented.

Partition function without monomers
which can be written for the special case of the square geometry where g L/2 is a number sequence (OEIS A065072) 19 [41]. Another observation is that the number of configuration on the square L × L is always even. It is less trivial to notice that {g p } is a sequence of odd numbers satisfying the relation [76] g p = p + 1(mod 32) if p even The exact solution of the dimer model with one boundary monomer allows for the same kind of number theory analysis (cf. [96] for details). The aim of the following sections is to look in more details at the form of the partition function of a dimer model of on an L × L square (L even) lattice with 2n monomers. One allows the 2n monomers to be anywhere on the four boundaries of the square (see Fig. 20).

Partition function with two boundary monomers
We saw previously that the expression of this partition function Q 2n is related to the pure dimer model Q 0 by the formula where the size of the matrix C depends on the number of monomers. Previously, the partition Q 0 has been shown to possess a remarkable expression (84) and we would like to determine whether or not, the partition function Q 2n admit the same kind of properties. In the case of two monomers anywhere on the boundaries, we saw that W = M and then the expression (86) reduces to

Inline boundary monomers
Initially, we choose to restrict the monomers to live on the same boundary m i = m j = L with n i , n j ∈ [1, L] (cf. Fig. 20(a)). In that particular situation the matrix elements C ij take the following form (53)  Table 3 Correlation function C ij for a boundary monomer (m i = m j = L) fixed on the first site n j = 1 as function of the ordinate n i for several system sizes L (see Fig. 20(a)). The value of C ij where the two monomers are on the corner (1, 1), (2, 1) is always equal to 1/2, because it is equivalent to force a dimer to be on the corner and then split the number of configuration by two. Bottom line: Values of the sequence g L for L = 4..14.
C ij (L) where = 0 if n i ∈ Z 2p and n j ∈ Z 2p+1 or conversely.
In Table 3, we evaluate this expression using MATHEMATICA ® , restricting one monomer to be in n j = 1 and the second to be between 1 to L for several system sizes. One can observe that there is a curious relation between the expression C ij and the sequence g L/2 present in the partition function Q 0 , more precisely one can deduce a proportionality relation which appears to be valid for all system sizes L in the case of inline monomers.

General case
The general expression of the matrix elements of the correlations between boundary monomers, Eq. (73), valid in all the geometries of Fig. 20 can be written as In Table 4, we evaluate this expression for the two other geometries. The same relation holds in this case as well, and we conjecture that the expression of the 2-point correlation takes the following form no matter the positions of the two monomers on the boundaries, where α (2) ij depends only on the positions of the two monomers and on the system size L. Consequently, the partition function of the dimer model with two boundary monomers reads Q 2 (L) = α (2) ij (L).g L/2 (93) Table 4 Top: Correlation function C ij between a monomer at position (m i = n i = 1) and another at position (m j = L, n j ) for several system sizes L and for n j = 1, 3, 5, 7, 9, 11 (cf. Fig. 20(b)). We notice that the expression of the last line in each row is the half of the expression of the penultimate line. Bottom: Correlation function C ij between opposite side monomers for several system sizes L and for n i = n j = 1, 2, 3, 4 (see Fig. 20(c)).
C ij (L)

Partition function with 2n boundary monomers
It is worth looking at higher number of monomers to conjecture a more general form of the partition function. We have conjectured that the matrix elements of the correlation matrix are proportional to the sequence g L/2 , thus thanks to the general pfaffian solution with 2n monomers equation (86), we obtain the formulas where the pure partition function takes the form (84), therefore the partition functions are proportional to power of g L/2 Q 0 (L) = 2 L/2 .g 2 L/2 , Q 2 (L) = α (2) ij (L).g L/2 , Q 4 (L) = α (4) ij kl (L).g 0 L/2 , Q 6 (L) = α (6) ij klmn (L).g −1 L/2 , . . . (95) which can be generalized for 2n monomers at positions i 1 , i 2 , . . . , i 2n in such a way that a relation between Q 2p and Q 2q can be found (p, q 1), dropping all the indices but p and q for simplicity, we found ex hypothesi valid for 2p and 2q monomers anywhere on the boundaries of the square lattice. Finally all these numerical relations between dimer partition functions with and without boundary monomers are the consequence of Eq. (86) and Eq. (92), which are unfortunately no longer valid for bulk monomers.

Conclusion
In this work the classical dimer model was discussed in great details both in a fermionic and bosonic field theory formulation. The bosonic formulation of the dimer model is based on the so-called height mapping and it is well suited for phenomenological predictions about correlations between dimers and monomers in a Coulomb gas context. Then we presented a practical and complete fermionic solution of the 2d dimer model on the square lattice with an arbitrary number of monomers. Furthermore, the Tzeng-Wu solution of the dimer model with a boundary monomer was found to be included in our theory. Interpretations of finite size effects of the Tzeng-Wu solution in a CFT/Coulomb gas framework has been performed, and we showed that a careful examination of boundary conditions in the model allowed us to recover the central charge of the free fermion/free boson field theory. The exact expression of correlation functions between monomers has been written in terms of the product of two pfaffians, and we gave an explicit formula for boundary correlations valid for the four boundaries of the rectangle. This solution has been used to compute correlations for several configurations in order to extract bulk, surface and corner scaling dimensions for dimer and monomer operators. All these results were interpreted in the Coulomb gas formalism, and we showed that all the predictions of the CFT were in accordance with a c = 1 theory. Last but not least, the exact closed-form expression of correlations between boundary monomers has been extensively used to extract some combinatorial and numerical informations about the partition function of the model. Furthermore, an unexpected relation has been found between partition functions with and without boundary monomers, and has been demonstrated thanks to our pfaffian solution. Generally, this Grassmann method can also be used for studying more general correlation functions, thermodynamical quantities, or transport phenomena of monomers. Other types of lattices, such as hexagonal lattice and other boundary conditions, can also be considered, as well as more precise comparisons with CFT results about rectangle geometry [129]. The same analysis of corner contribution to free energy as well as critical exponents can be studied in the interacting dimer model using the height mapping and results will be presented elsewhere. A future challenge emerging out of this present work is the study of other two dimensional dimer related models as the trimer model [53] or the four-color model [91,46] which can be seen as an interacting colored dimer model. Work in those directions is in progress. mapping and Christophe Chatelain, Jérôme Dubail, Malte Henkel, Jesper Jacobsen, Shahin Rouhani and Loïc Turban for insightful discussions and useful comments at different stages of the elaboration of this article. This work was partly supported by the Collège Doctoral Leipzig-Nancy-Coventry-Lviv (Statistical Physics of Complex Systems) of UFA-DFH. which leads to the generalization It is obvious that this definition fulfills the constraint of translational invariance which requires Changes of coordinates are required to preserve the anticommuting structure of the Grassmann algebra, this allows non-singular linear transformations of the form b i = j A ij a i . One then can verify that by setting f (a) = F (b) one can obtain the following relation at variance with the commuting case in which the factor on the right hand side would have been det −1 A. We define D[a, ā] = i da i dā i the Grassmann measure. In the multidimensional integral, the symbols da 1 , . . . , da N are again anticommuting with each other. The basic expression of the Grassmann analysis concerns the Gaussian fermionic integrals [127] which is related to the determinant where {a i , ā i } is a set of completely anticommuting Grassmann variables, the matrix in the exponential is arbitrary. The two Grassmann variables a i and ā i are independent and not conjugate to each other, they can been seen as component of a complex Grassmann variables. The Gaussian integral of the second kind is related to the pfaffian of the associated skew-symmetric matrix The pfaffian form is a combinatorial polynomial in A ij , known in mathematics for a long time.
The pfaffian and determinant of the associated skew-symmetric matrix are algebraically related by det A = (pfA) 2 . This relation can be most easily proved in terms of the fermionic integrals. The linear superpositions of Grassmann variables are still Grassmann variables and it is possible to make a linear change of variables in the integrals. The only difference with the rules of the common analysis, is that the Jacobian will now appear in the inverse power. New variables of integration can be introduced, in particular, by means of the transformation to the momentum space. The permanent of A and the so-called haffnian can be written with Grassmann variables as well which are connected by the formula permA = (hfA) 2 . We recall that the definition of the permanent differs from that of the determinant in that the signatures of the permutations are not taken into account.

Appendix B. Plechko mirror symmetry
In this appendix we briefly recall the method of resolution of the 2d dimer model based on the integration over Grassmann variables and factorization principles for the partition function introduced in the context of 2d Ising model [117]. The general partition function for a graph with N vertices can be written, for a square lattice of size L × L with L even, as L m,n (1 + t x η mn η m+1n )(1 + t y η mn η mn+1 ), where η mn are nilpotent and commuting variables on every vertices of the square lattice. The integrals can be done if we introduce a set of Grassmann variables (a mn , ā mn , b mn , b mn ) (cf. Fig. 21(a)), such that This decomposition allows for an integration over variables η mn , after rearranging the different link variables A mn := 1 + a mn η mn , Ā m+1n := 1 + t xāmn η m+1n , B mn := 1 + b mn η mn and B mn+1 := 1 + t ybmn η mn+1 . Then the partition function becomes where we use the notation for the measure of integration Tr {a,ā,b,b,η} X (a,ā, b,
Then, the non-commuting link variables are moved in such a way that each η mn is isolated and can be integrated directly. This rearrangement is possible in two dimensions thanks to the mirror ordering introduced by Plechko for the Ising model. The ordering process can be detailed as follows where the products are ordered according to the orientation of the arrows. The Grassmann terms in brackets (· · ·) on the first line of the previous equation are commuting objects, since they are integral representations of commuting scalars. This also imposes the boundary conditions Ā 1n = 1, Ā L+1n = 1, B m1 = 1, and B mL+1 = 1, or ā 0n =ā Ln =b m0 =b mL = 0 (for open boundary conditions only). We finally obtain the following exact expression The integration over the η mn variables is performed exactly, recursively from m = 1 to m = L for each n. Each integration leads to a quantity L mn = a mn + b mn + t xām−1n + (−1) m+1 t ybmn−1 which is moved to the left of the products over m, hence a minus sign is needed in front of b each time an L mn crosses the product of B terms on the left. Finally becomes an integration over products of linear Grassmann terms. This can be further simplified by introducing additional Grassmann variables c mn such that L mn = dc mn exp(c mn L mn ). (118) This expresses Q 0 as a Gaussian integral over variables (a, ā, b, b , c), and therefore Q 0 is a simple determinant of a quadratic form. Indeed, after partially integrating over variables (a, ā, b, b ) and symmetrization of the expressions, one obtains with c 0n = c L+1n = c m0 = c mL+1 = 0. Inserting Eq. (120) into Eq. (119), and using two following sum identities where, for example, c −pq is a short notation for c L+1−pq . The summation is performed only for 1/4 of the Fourier modes (see Fig. 21(b)), since the other are related by the symmetry p → L + 1 − p and q → L + 1 − q. This block representation is convenient for computing the remaining integrals over the momenta, as a product of cosine functions as found by Kasteleyn, Temperley and Fischer 4t 2 x cos 2 πp L + 1 + 4t 2 y cos 2 πq L + 1 . The main ingredient of this recipe is the mirror factorization equation (116) of the partition function, which allows a direct integration over nilpotent variables. This mirror factorization is then the major obstacle to the generalization of this formalism in d > 2. Indeed, for the 3d case, a generalization of this factorization is far from obvious and remains to find.

Appendix C. General monomer-dimer partition function
The general monomer-dimer problem is a much more complex and challenging problem in statistical physics and combinatorics, because the position of the monomers are not fixed either than their number (cf. Fig. 22 for L = 2). From the point of view of theoretical physics, the number of monomers divided by the number of occupied site defines the monomer density ρ. It is long known that the full phase diagram of the monomer-dimer model does not admit any phase transition for ρ > 0 [59]. Furthermore the behavior of monomer-monomer correlations for finite density has been studied numerically [101], and strong evidences for exponential correlations have been established, in accordance with mean-field calculations using Grassmann variables [114]. From a computational point of view, the problem has been shown to belong to the #P -complete enumeration class [75] and all the methods available are either efficient but approximative [8,82] or exact but desperately slow [1]. In this short appendix one shows how to use our exact solution to express the partition function of this enumerative problem.
This number grows as 2 L 2 when the size of the lattice goes to infinity, making the problem impossible to solve analytically. Since our method allows to calculate exactly the partition function of the dimer model on a square lattice of size M × N with an arbitrary even number of monomers, then we can formally write down the full monomer-dimer partition function as a sum over the number and the positions of monomers which becomes simpler in the boundary case