We review free fermion, melting crystal, and matrix model representations of wall-crossing phenomena on local, toric Calabi-Yau manifolds. We consider both unrefined and refined BPS counting of closed BPS states involving D2- and D0-branes bound to a D6-brane, as well as open BPS states involving open D2-branes ending on an additional D4-brane. Appropriate limit of these constructions provides, among the others, matrix model representation of refined and unrefined topological string amplitudes.

1. Introduction

This paper is devoted to some aspects of counting of BPS states in a system of D𝑝-branes, with even 𝑝, in type IIA string compactifications. The problems of BPS counting span a vast area of research in supersymmetric gauge and string theories. Their important feature is a special, nonconstant character of BPS multiplicities: their values depend on various moduli and jump discontinuously along some special loci in the corresponding moduli space, so called walls of marginal stability. The pattern of these jumps follows wall-crossing formulas, found from physical perspective by Denef and Moore [1] and, in more general context, formulated mathematically by Kontsevich and Soibelman [2]. The regions of the moduli space in between walls of marginal stability, in which BPS multiplicities are (locally) constant, are called chambers.

The BPS states we are interested in, and which we will refer to as closed BPS states, arise as bound states of a single D6-brane with arbitrary number of D0 and D2-branes wrapping cycles of a toric Calabi-Yau space. More generally, we will also consider open BPS states, which arise when an additional D4-brane spans a Lagrangian submanifold inside the Calabi-Yau space and supports open D2-branes attached to it. The closed and open BPS states give rise, respectively, to single-particle states in the effective four-dimensional and two-dimensional theory (in remaining, space-time filling directions of, resp., D6 and D4-branes). In this context, the character of BPS multiplicities can be understood in much detail, and it relates to other interesting exactly solvable models: free fermions, crystal, and matrix models. In brief, these connections arise as follows. Firstly, BPS states we consider turn out to be in one-to-one correspondence with configurations of certain statistical models of melting crystals. The structure of these crystals depends on geometry of the underlying Calabi-Yau space, as well as on the chamber one is considering. In consequence, BPS counting functions, upon appropriate identification of parameters, coincide with generating functions of melting crystals. It turns out that the structure of these crystals can be given a free fermion representation. Furthermore, once such free fermion formulation is known, it can also be represented in terms of matrix models. Connection with vast theory of matrix models has many interesting mathematical and physical consequences and allows to shed new light on wall-crossing phenomena. The aim of this paper is to explain these connections.

The BPS generating functions which we consider are intimately related to topological string amplitudes on corresponding Calabi-Yau spaces. This relation is most transparent in the physical derivation discussed in Section 2, which relies on lifting the D-brane system to M-theory. The M-theory viewpoint makes contact with original formulation of closed topological strings by Gopakumar and Vafa [3, 4], and open topological strings by Ooguri and Vafa [5]. In particular, in one specific, so-called noncommutative chamber, the BPS-generating function is given as the modulus square of the topological string partition function. In all other chambers, BPS generating functions can be uniquely determined from that noncommutative result. There is also another special, so-called commutative chamber, in which BPS generating function coincides (up to the factor of MacMahon function) with the topological string partition function. For toric manifolds which we consider, such topological string amplitudes can be constructed, among the other, by means of the powerful topological vertex formalism [6]. Relation to crystal models was in fact first understood in this topological string chamber [7–9]. One advantage of the formalism presented in this paper is the fact that it allows to construct matrix model representation of all these generating functions (so, in particular, matrix model representation of topological string amplitudes).

In more detail, we will consider generating functions of D2 and D0-branes bound to a single D6-brane of the following form: 𝒡BPSξ€·π‘žπ‘ ξ€Έ=,𝑄𝛼,𝛽Ω(𝛼,𝛽)π‘žπ›Όπ‘ π‘„π›½,(1.1) where π›Όβˆˆβ„€ is D0-brane charge and π›½βˆˆπ»2(𝑋,β„€) is D2-brane charge. Multiplicities Ξ©(𝛼,𝛽) jump when central charges (which itself are functions of KΓ€hler moduli) of building blocks of a bound state align, and therefore these generating functions are locally constant functions of KΓ€hler moduli. Along the walls of marginal stability, the degeneracies Ξ©(𝛼,𝛽) change and indeed obey wall-crossing formulas of [1, 2] mentioned above.

If there is an additional D4-brane which spans a Lagrangian submanifold inside the Calabi-Yau space, in addition to the above closed BPS states, one can consider also open BPS states of D2-branes with boundaries ending on a one-cycle 𝛾 on this D4-brane. In this case, the BPS states arise on the remaining two-dimensional world-volume of the D4-brane. The holonomy of the gauge field along 𝛾 provides another generating parameter 𝑧, so that open BPS-generating functions take form 𝒡openBPSξ€·π‘žπ‘ ξ€Έ=,𝑄𝛼,𝛽,𝛾Ω(𝛼,𝛽)π‘žπ›Όπ‘ π‘„π›½π‘§π›Ύ.(1.2) As we will show, generating functions of such open BPS states can be identified with integrands of matrix models mentioned above.

One more important aspect of BPS counting is referred to as refinement, and amounts to refining BPS counting by introducing one more parameter, customarily denoted 𝛽. The refinement can be introduced from several perspectives which give rise to identical results; however, their fundamental common origin is still not fully understood. We will introduce refinement by distinct counting of states with different π‘†π‘ˆ(2) spins inside spacetime 𝑆𝑂(4) rotation group in the generating function (1.1). In [10], it was argued that this physical viewpoint should agree with the mathematical counterpart of motivic deformation [2], and also a refined version of a crystal model was constructed. Another notion of refinement arises in Nekrasov partition functions, which are defined in a nontrivial gravitational (so-called Ξ©-) background parametrized by two parameters πœ–1 and πœ–2 [11]. Nekrasov partition functions can also be defined for five-dimensional gauge theories and then they agree with topological string amplitudes. In particular, the formalism of the topological vertex [6] has also been extended to the refined context in [12], and shown to reproduce relevant Nekrasov partition functions. Also BPS generating functions, in the limit of commutative chamber, are known to reproduce refined topological string amplitudes with 𝛽=βˆ’πœ–1/πœ–2 [13]. However, the worldsheet definition of refined topological string amplitudes is not fully understood.

As an exemplary and, hopefully, inspiring application of the entire formalism presented in this paper, in the final Section 6 we derive matrix model representation of the refined topological string partition function for the conifold. The refined matrix model which we find has a standard measure; however, its potential is deformed by 𝛽-dependent terms. It is obtained by constructing appropriate refined crystal model and free fermion representation, and subsequently reformulating this representation in matrix model form. Finally, taking the limit of the commutative chamber, we obtain matrix model representation of the refined topological string amplitude. Even though we demonstrate this result in the conifold case, with some technical effort it can be generalized to other toric manifolds which we consider (As we recall in Section 6, refined topological string amplitudes were also postulated to be reproduced by another type of matrix models, so-called 𝛽-deformed ones (whose Vandermonde measure is deformed by raising it to power 𝛽); however, explicit computations showed that this cannot be the correct representation of refined amplitudes.).

1.1. Short Literature Guide

The literature on the topics presented in this paper is extensive and still growing, and we unavoidably mention just a fraction of important developments. The relation between Donaldson-Thomas invariants for the noncommutative chamber of the conifold was first found by SzendrΕ‘i [14]. It was generalized to orbifolds of β„‚3, and related to free fermion formalism, by Bryan and Young [15]. The relation to free fermions and crystals was extended to a large class of toric manifolds without compact four-cycles [16, 17]. These developments were accompanied by other mathematical works [18, 19].

In parallel to the above-mentioned mathematical activity, wall-crossing phenomena for local Calabi-Yau manifolds were analyzed from physical viewpoint. The analysis of nontrivial BPS counting for the conifold was described by Jafferis and Moore in [20]. This and more general cases were related to quivers and crystal models in [21, 22]. Derivation of BPS degeneracies from M-theory viewpoint and relation to closed topological strings were discussed in [23], and generalized to open BPS counting in [24–27]. Relations to matrix models, discussed for plane partitions with some other motivation in [28], were extended to other crystal models relevant for BPS counting in [29], and also in [30]. Subsequently, it was related to open BPS counting in [27]. Refined BPS counting was related to crystal models in [10, 13], and corresponding matrix models were constructed in [31].

Let us also mention some other, related works devoted to crystals and free fermions. The fermionic construction of MacMahon function for β„‚3 was originally presented in [7], and its relation to open topological strings and more complicated Calabi-Yau manifolds were discussed in [32–34]. Newer ideas, analyzing more complicated systems involving D4-branes, were presented in [35, 36]. More expository presentations of various aspects described here can be found in [37, 38]. A general introduction to mathematical and physical aspects of mirror symmetry can be found in [39].

1.2. Plan

The plan of this paper is as follows. In Section 2, we introduce BPS generating functions and present one possible derivation of their form, which relies on the M-theory interpretation of a D-brane system, following [23–25, 27]. In Section 3, we provide a little mathematical background and introduce notation pertaining to toric Calabi-Yau manifolds, free fermion formalism, and matrix models. In Section 4, we introduce fermionic formalism for BPS generating functions and present corresponding crystal models, building on earlier ideas of [7, 15] and following [16]. In Section 5, we reformulate the problem of closed BPS counting in terms of matrix models and relate it to open BPS counting [27, 29]. In Section 6, we refine our analysis, present refined BPS generating functions and crystals [10], and construct corresponding refined matrix models [31].

2. BPS Generating Functions

In this section, we introduce generating functions of BPS states of D-branes in toric Calabi-Yau manifolds. Our task in the rest of this paper is to provide interpretation of these generating functions in terms of free fermions, melting crystals, and matrix models. These generating functions can be derived using wall-crossing formulas, as was done first in the unrefined [20] and refined [10] conifold case, and later generalized to arbitrary geometry without compact four-cycles in [18, 19]. On the other hand, we will focus on a simpler physical derivation of BPS generating functions which uses the lift of the D-brane system to M-theory [23]. This also makes contact with M-theory interpretation of topological string theory and allows to express BPS counting functions in terms of topological string amplitudes. Moreover, this M-theory derivation can be extended to the counting of open BPS states,that is, open D2-branes attached to additional D4-brane, which we are also interested in [24, 25, 27].

We start this section by reviewing the M-theory derivation of (unrefined) closed and open BPS generating functions. Then, to get acquainted with a crystal interpretation of these generating functions, we discuss their crystal interpretation in simple cases of β„‚3 and conifold. Later, using fermionic interpretation, we will generalize this crystal representation to a large class of toric manifolds without compact four cycles.

2.1. M-Theory Derivation

We start by considering a system of D2 and D0-branes bound to a single D6-brane in type IIA string theory. It can be reinterpreted in M-theory as follows [23]. When additional 𝑆1 is introduced as the eleventh dimension transversely to the D6-brane, then this D6-brane transforms into a geometric background of a Taub-NUT space with unit charge [40]. The Taub-NUT space is a circle fibration over ℝ3, with a circle 𝑆1TN attaining a fixed radius 𝑅 at infinity, and shrinking to a point in the location of the original D6-brane. From M-theory perspective, bound states involving D2 and D0-branes are interpreted as M2-branes with momentum on a circle. Therefore, the counting of original bound states to the D6-brane is reinterpreted as the counting of BPS states of M2-branes in the Taub-NUT space. While in general this is still a nontrivial problem, for the purpose of counting BPS degeneracies we can take advantage of their invariance under continuous deformations of the Taub-NUT space, in particular under deformations of the radius 𝑅. We can therefore consider taking this radius to infinity, whereupon BPS counting is reinterpreted in terms of a gas of particles in ℝ5. To make the problem fully tractable, we have to ensure that the particles are noninteracting, which would be the case if moduli of the Calabi-Yau would be tuned so that M2-branes wrapped in various ways would have aligned central charges. This can be achieved when KΓ€hler parameters of the Calabi-Yau space are tuned to zero. However, to avoid generation of massless states, at the same time one has to include nontrivial fluxes of the M-theory three-form field through the two cycles of the Calabi-Yau and 𝑆1TN. In type IIA, this results in the 𝐡-field flux 𝐡 through two cycles of Calabi-Yau. Finally, to avoid creation of the string states arising from M5-branes wrapping four cycles in Calabi-Yau, we simply restrict considerations to manifolds without compact four cycles. For a state arising from D2-brane wrapping a class 𝛽, the central charge then reads 1𝑍(𝑙,𝛽)=𝑅(𝑙+𝐡⋅𝛽),(2.1) where 𝑙 counts the D0-brane charge, which is taken positive to preserve the same supersymmetry.

Under the above conditions, the counting of D6-D2-D0 bound states is reinterpreted in terms of a gas of particles arising from M2-branes wrapped on cycles 𝛽. The excitations of these particles in ℝ4, parametrized by two complex variables 𝑧1,𝑧2, are accounted for by the modes of the holomorphic field Φ𝑧1,𝑧2ξ€Έ=𝑙1,𝑙2𝛼𝑙1,𝑙2𝑧𝑙11𝑧𝑙22.(2.2) Decomposing the isometry group of ℝ4 as 𝑆𝑂(4)=π‘†π‘ˆ(2)Γ—π‘†π‘ˆ(2)β€², there are π‘π›½π‘š,π‘šξ…ž five-dimensional BPS states of intrinsic spin (π‘š,π‘šβ€²). We are interested in their net number arising from tracing over π‘†π‘ˆ(2)β€² spins π‘π‘šπ›½=ξ“π‘šβ€²(βˆ’1)π‘šβ€²π‘π‘š,π‘šβ€²π›½.(2.3) The total angular momentum of a given state contributing to the index is 𝑙=𝑙1+𝑙2+π‘š. Finally, in a chamber specified by the moduli 𝑅 and 𝐡, the invariant degeneracies can be expressed as the trace over the corresponding Fock space 𝒡BPS=ξ‚€TrFockπ‘žπ‘„0𝑠𝑄𝑄2|||chamber=𝛽,π‘šξ‘π‘™1+𝑙2=𝑙1βˆ’π‘žπ‘™1+𝑙2𝑠+π‘šπ‘„π›½ξ‚π‘π‘šπ›½||||chamber=ξ‘βˆžπ›½,π‘šξ‘π‘™=1ξ€·1βˆ’π‘žπ‘ π‘™+π‘šπ‘„π›½ξ€Έπ‘™π‘π‘šπ›½|||chamber,(2.4) where the subscript chamber denotes restriction to those factors in the above product, which represent states which are mutually BPS 𝑍(𝑙,𝛽)>0βŸΊπ‘žπ‘ π‘™+π‘šπ‘„π›½<1.(2.5) As usual, 𝑄=π‘’βˆ’π‘‡ and π‘žπ‘ =π‘’βˆ’π‘”π‘  above encode, respectively, the KΓ€hler class 𝑇 and the string coupling 𝑔𝑠 (we wish to distinguish carefully π‘žπ‘  which encodes string coupling, from a counting parameter π‘ž which will arise in what follows in crystal interpretation). The above condition on central charges is crucial in determining a particular form of the BPS generating functions. If we would restrict products in the formula (2.4) to factors with only positive 𝛽, we would get (up to possibly some factor of MacMahon function) the Gopakumar-Vafa representation of the topological string amplitude. With all negative and positive values of 𝛽, we would get modulus square of the topological string partition function. Therefore, the upshot of [23] is that in general the above BPS generating function can be expressed in terms of the closed topological string partition function 𝒡BPS=𝒡top(𝑄)𝒡topξ€·π‘„βˆ’1ξ€Έ||chamber,(2.6) where chamber restriction is to be understood as picking up only those factors in Gopakumar-Vafa product representation of 𝒡top for which (2.5) is satisfied. In this context, we will often refer to the choice of a chamber as a closed BPS chamber. The (instanton part of the) closed topological string partition function entering the above expression is given by [3, 4] 𝒡topξ€·π‘ž(𝑄)=π‘€π‘ ξ€Έβˆžπœ’/2𝑙=1𝛽>0,π‘šξ€·1βˆ’π‘„π›½π‘žπ‘ π‘š+π‘™ξ€Έπ‘™π‘π‘šπ›½,(2.7) where 𝑀(π‘žπ‘ βˆ)=𝑙(1βˆ’π‘žπ‘™π‘ )βˆ’π‘™ is the MacMahon function and πœ’ is the Euler characteristic of the Calabi-Yau manifold.

To be more precise, an identification as a topological string partition function or its square arises if 𝑅>0 in (2.1). Because 𝑅 arises just as a multiplicative factor in (2.1), degeneracies depend only on its sign. Therefore, another extreme case corresponds to negative 𝑅 and 𝐡 sufficiently small, when only a single D6-brane contributes to the partition function 𝒡(𝑅<0,0<𝐡β‰ͺ1)=1.(2.8) More generally, for 𝑅<0, BPS generating functions often (but not always) take finite form.

In what follows we denote BPS generating functions in chambers with positive 𝑅 by 𝒡, and in chambers with negative 𝑅 by 𝒡 (and often omit the subscript BPS). Topological string partition functions will be denoted by 𝒡top, while generating functions of melting crystals by ordinary 𝑍.

The above structure can be generalized by including in the initial D6-D2-D0 configuration additional D4-branes wrapping Lagrangian cycles in the internal Calabi-Yau manifold and extending in two space-time dimensions [24, 25, 27]. For simplicity, we consider a system with a single D4-brane wrapping a Lagrangian cycle. There are now additional BPS states in two remaining spacetime dimensions arising from open D2-branes ending on these D4-branes. Their net degeneracies 𝑁𝑠,𝛽,𝛾 are characterized, firstly, by the 𝑆𝑂(2) spin 𝑠 whose origin is most clearly seen from the M-theory perspective [5, 41]. Secondly, they depend on two-cycles 𝛽 wrapped by open M2-branes, as well as one-cycles 𝛾 on which these M2-branes can end (In case of 𝑁 D4-branes wrapping the same Lagrangian cycle, these states would additionally arise in representations 𝑅 of π‘ˆ(𝑁) [5]. In case of a single brane, this reduces to π‘ˆ(1), and such a dependence can be reabsorbed into a parameter specifying a choice of 𝛾.).

Lifting this system to M-theory, we obtain a background of the form Taub-NUTΓ—Calabi-Yau×𝑆1, with the additional D4-brane promoted to M5-brane. This M5-brane wraps the Lagrangian submanifold 𝐿 inside Calabi-Yau, the time circle 𝑆1, and ℝ+×𝑆1TN inside the Taub-NUT space. A part of this Lagrangian 𝐿 is a torus 𝑇2=𝑆1TN×𝑆1, which will lead to some modular properties of the BPS counting functions: this modularity will be manifest in one chamber, where the open topological string amplitude will be completed to the product of πœƒ functions. This M5-brane also breaks the 𝑆𝑂(4) spatial symmetry down to 𝑆𝑂(2)×𝑆𝑂(2)ξ…ž. We denote the spins associated to both 𝑆𝑂(2) factors, respectively, by 𝜎 and πœŽβ€², and the degeneracies of particles with such spins by π‘πœŽ,πœŽξ…žπ›½,𝛾. In addition to closed KΓ€hler parameters 𝑄=π‘’βˆ’π‘‡, let us also introduce open ones related to discs wrapped by M2-branes 𝑧=π‘’βˆ’π‘‘. The real and imaginary parts of 𝑇 encode, respectively, the sizes of two-cycles 𝛽 and the value of the 𝐡-field through them. The real and imaginary parts of 𝑑 encode, respectively, sizes of the discs and holonomies of the gauge fields around them. Similarly as in the closed string case, to get nontrivial ensemble of mutually supersymmetric states, we set the real parts of 𝑇 and 𝑑 to zero, and consider nontrivial imaginary parts.

From the M-theory perspective, we are interested in counting the net degeneracies of M2-branes ending on this M5-brane π‘πœŽ,𝛽,𝛾=ξ“πœŽβ€²(βˆ’1)πœŽβ€²π‘πœŽ,πœŽβ€²π›½,𝛾.(2.9) In the remaining three-dimensional space, in the π‘…β†’βˆž limit, the M2-branes ending on the M5-brane are represented by a gas of free particles. These particles have excitations in ℝ2 which we identify with the 𝑧1-plane. To each such BPS particle, similarly as in the closed string case discussed above and in [23, 40], we can associate a holomorphic field Φ𝑧1ξ€Έ=𝑙𝛼𝑙𝑧𝑙1.(2.10) The modes of this field create states with the intrinsic spin 𝑠 and the orbital momentum 𝑙 in the ℝ2 plane. The derivation of the BPS degeneracies relies on the identification of this total momentum 𝜎+𝑙 in the π‘…β†’βˆž limit, with the Kaluza-Klein modes associated to the rotations along 𝑆1TN for the finite 𝑅, following the five-dimensional discussion in [40, 42].

The BPS generating functions we are after are given by a trace over the Fock space built by the oscillators of the second quantized field Ξ¦(𝑧1) and restricted to the states which are mutually supersymmetric. In such a trace, each oscillator from (2.10) gives rise to one factor of the form (1βˆ’π‘žπ‘ πœŽ+π‘™βˆ’1/2𝑄𝛽𝑧𝛾)Β±1, where the exponent Β±1 corresponds to the bosonic or fermionic character of the top component of the BPS state, 𝒡openBPS=ξ‘βˆžπœŽ,𝛽,𝛾𝑙=1ξ€·1βˆ’π‘žπ‘ πœŽ+π‘™βˆ’1/2π‘„π›½π‘§π›Ύξ€Έπ‘πœŽ,𝛽,𝛾|||chamber,(2.11) where the product is over either both positive or both negative (𝛽,𝛾). The parametersπ‘ž,𝑄 and 𝑧 specify the chamber structure: the restriction to a given chamber is implemented by imposing the condition on a central charge, analogous to (2.5), π‘žπ‘ πœŽ+π‘™βˆ’1/2𝑄𝛽𝑧𝛾<1.(2.12) This condition in fact specifies a choice of both closed and open chambers. The walls of marginal stability between chambers correspond to subspaces where, for some oscillator, the above product becomes 1, and then the contribution from such an oscillator drops out from the BPS generating function.

Similarly as in the closed string case, the above degeneracies can be related to open topological string amplitudes, rewritten in [5] in the form 𝒡opentop=expβˆžξ“π‘›=1ξ“πœŽξ“π›½,𝛾>0π‘πœŽ,𝛽,π›Ύπ‘žπ‘ π‘›πœŽπ‘„π‘›π›½π‘§π‘›π›Ύπ‘›ξ€·π‘žπ‘ π‘›/2βˆ’π‘žπ‘ βˆ’π‘›/2ξ€Έξƒͺ,(2.13) with integer Ooguri-Vafa invariants π‘πœŽ,𝛽,𝛾 (In case of 𝑁 D4-branes wrapping a Lagrangian cycle, this structure is again more complicated, because the states in ℝ3 arise in representations of π‘ˆ(𝑁) [5]. This requires replacing the factor 𝑧𝑛𝛾 by the sum βˆ‘π‘…Tr𝑅𝑉𝑛 of traces in all possible representations 𝑅 of this π‘ˆ(𝑁) of the matrix 𝑉 encoding holonomies of the gauge fields. For simplicity we restrict here to the simplest case.). This formula represents in fact a series of quantum dilogarithms 𝐿𝑧,π‘žπ‘ ξ€Έξƒ©ξ“=exp𝑛>0π‘§π‘›π‘›ξ€·π‘žπ‘ π‘›/2βˆ’π‘žπ‘ βˆ’π‘›/2ξ€Έξƒͺ=βˆžξ‘π‘›=1ξ€·1βˆ’π‘§π‘žπ‘ π‘›βˆ’1/2ξ€Έ(2.14) and can be written in the product form 𝒡opentop(𝑄,𝑧)=πœŽξ‘βˆžπ›½,𝛾>0𝑛=1ξ€·1βˆ’π‘„π›½π‘§π›Ύπ‘žπ‘ πœŽ+π‘›βˆ’1/2ξ€Έπ‘πœŽ,𝛽,𝛾.(2.15) Comparing with (2.11) we conclude that the BPS counting functions take form of the modulus square of the open topological string amplitude 𝒡openBPS=𝒡opentop(𝑄,𝑧)𝒡opentopξ€·π‘„βˆ’1,π‘§βˆ’1ξ€Έ||chamber.(2.16)

Similarly as in the closed string case, there are also a few particularly interesting chambers to consider. For example, in the extreme chamber corresponding to Im 𝑇, Im 𝑑→0, the trace is performed over the full Fock space and yields the modulus square of the open topological string partition function. In this case, the quantum dilogarithms arise in pairs, which (using the Jacobi triple product identity) combine to the modular function πœƒ3/πœ‚; in consequence, the total BPS generating function is modular and expressed as a product of such functions.

2.2. Crystal Interpretation

Closed BPS generating functions (2.4) turn out to be generating functions of statistical models of crystals, when parameters relevant for both interpretations are appropriately matched. Physical reasons for such relations have been given in [8, 21, 22], and mathematical interpretation arose from works [9, 14, 15]. Such crystal interpretation arises also from the fermionic formulation [16, 17], as we will review below. These crystals, in a more intricate way [27], encode also open BPS generating functions (2.11). However, before discussing details of all these constructions, in this introductory section we present crystal models for two simplest toric Calabi-Yau manifolds, that is, β„‚3 and conifold.

β„‚3 is the simplest Calabi-Yau manifold. It has no compact two-cycles, so relevant BPS states are bound states of arbitrary number of D0-branes with a single D6-brane wrapping entire β„‚3. Their generating function is therefore expressed in terms of a single parameter π‘žπ‘ =π‘’βˆ’π‘”π‘ . There is just a single nonzero Gopakumar-Vafa invariant 𝑁0𝛽=0=βˆ’1, and as follows from (2.4) this generating function coincides with the so-called MacMahon function 𝒡BPS=βˆžξ‘π‘™=11ξ€·1βˆ’π‘žπ‘™π‘ ξ€Έπ‘™ξ€·π‘ž=𝑀𝑠.(2.17)

On the other hand, the MacMahon function is a generating function of plane partitions, that is, three-dimensional generalization of Young diagrams. These plane partitions represent the simplest three-dimensional crystal model, namely, they can be identified with stacks of unit cubes filling the positive octant of ℝ3 space, as shown in Figure 1. A unit cube located in position (𝐼,𝐽,𝐾) can evaporate from this crystal only if all other cubes with coordinates (𝑖≀𝐼,𝑗≀𝐽,π‘˜β‰€πΎ) are already missing. A plane partition πœ‹ is weighted by the number of boxes it consists of |πœ‹|, with a weight π‘ž associated to a single box, so indeed 𝑍=πœ‹π‘ž|πœ‹|=βˆžξ“π‘™=0𝑝(𝑙)π‘žπ‘™=1+π‘ž+3π‘ž2+6π‘ž3+13π‘ž4+β‹―=𝑀(π‘ž),(2.18) where 𝑝(𝑙) is the number of plane partitions which consist of 𝑙 cubes. Therefore, plane partition generating function coincides with the BPS counting function 𝑍=𝒡BPS when a simple identification π‘žπ‘ =π‘ž(2.19) is made. From (2.6), it follows that the topological string partition function for β„‚3 is given by the square root of the MacMahon function 𝒡topξ€·π‘ž=𝑀𝑠1/2,(2.20) which is indeed true. The relevance of the MacMahon function for β„‚3 geometry was noticed for the first time in [3], and a statistical model interpretation of this result was proposed in [7].

The conifold provides another simple, yet nontrivial example of toric Calabi-Yau manifold. It consists of two β„‚3 patches glued intπ’ͺ(βˆ’1)βŠ•π’ͺ(βˆ’1)β†’β„™1, and it has one KΓ€hler class representing β„™1, parametrized by 𝑄=π‘’βˆ’π‘‡. This class can be wrapped by D2-branes, which bind with D0-branes to an underlying D6-brane and give rise to BPS states in low energy theory. In this case, there is already a nontrivial structure of chambers and walls, which was analyzed in [14, 18, 20, 21]. This structure is consistent with M-theory derivation discussed in Section 2.1. The generating functions of D6-D2-D0 bound states are parametrized by 𝑄 and π‘žπ‘ , and therefore corresponding crystal models consist of two-colored three-dimensional partitions. The KΓ€hler moduli space consists of several infinite countable sets of chambers, and in each chamber relevant crystal configurations take form of so-called pyramid partitions. These partitions are infinite or finite (resp. for positive and negative 𝑅 in (2.1)) and their size depends on the value of the 𝐡-field. This size changes discretely and the pyramid is enlarged when the value of the 𝐡-field crosses integer numbers, which changes the chamber in the moduli space, as explained in Section 2.1. Examples of such infinite pyramid partitions are given in Figure 2, and finite ones in Figure 3.

To write down explicitly BPS generating functions for the conifold in various chambers, we can take advantage of their relation to the topological string amplitude (2.6). The topological string partition function in this case reads 𝒡conifoldtopξ€·π‘ž(𝑄)=π‘€π‘ ξ€Έξ‘π‘˜β‰₯1ξ€·1βˆ’π‘„π‘žπ‘˜π‘ ξ€Έπ‘˜,(2.21) with the MacMahon function defined in (2.17). From this topological string partition function we can read off Gopakumar-Vafa invariants [3, 4] 𝑁0𝛽=0=βˆ’2,𝑁0𝛽=Β±1=1.(2.22)

Using the relation (2.6), we can now present conifold closed BPS generating functions in several sets of chambers. In the first set of chambers, we consider 𝑅>0 and positive 𝐡∈]𝑛,𝑛+1[ (for 𝑛β‰₯0). Firstly, for small 𝐡, there is so-called noncommutative chamber discussed first by SzendrΕ‘i [14], which corresponds to 𝑛=0. In this case, the pyramid crystal has just a single ball in the top row, as in Figure 2(a), and the BPS generating function is given by the square of the topological amplitude. On the other hand, for large 𝐡, that is π‘›β†’βˆž, we reach commutative chamber in which the length of the top row extends to infinity. In this case, the BPS generating function agrees, up to a single factor of MacMahon function, with the topological string amplitude. In between, there are chambers with 𝑛+1 balls in the top row, for which 𝒡conifoldπ‘›ξ€·π‘ž=𝑀𝑠2ξ‘π‘˜β‰₯1ξ€·1βˆ’π‘„π‘žπ‘˜π‘ ξ€Έπ‘˜ξ‘π‘˜β‰₯𝑛+1ξ€·1βˆ’π‘„βˆ’1π‘žπ‘˜π‘ ξ€Έπ‘˜.(2.23) These BPS generating functions are related to pyramid generating functions with two colors π‘ž0 and π‘ž1 upon the identification (which generalizes (2.19) in β„‚3 case) 𝒡conifold𝑛chambersβˆΆπ‘žπ‘ =π‘ž0π‘ž1,𝑄=βˆ’π‘žπ‘›π‘ π‘ž1.(2.24) Indeed, with this identification, the above counting functions agree with those of two-colored pyramid crystals with 𝑛+1 yellow balls in its top row 𝑍pyramidπ‘›ξ€·π‘ž0,π‘ž1ξ€Έξ€·π‘ž=𝑀0π‘ž1ξ€Έ2ξ‘π‘˜β‰₯𝑛+1ξ€·1+π‘žπ‘˜0π‘ž1π‘˜+1ξ€Έπ‘˜βˆ’π‘›ξ‘π‘˜β‰₯1ξ€·1+π‘žπ‘˜0π‘ž1π‘˜βˆ’1ξ€Έπ‘˜+𝑛.(2.25)

In the second set of chambers, we have 𝑅<0 and positive 𝐡∈]π‘›βˆ’1,𝑛[ (for 𝑛β‰₯1). It extends between the core region with a single D6-brane (2.8) and the chamber characterized by so-called Pandharipande-Thomas invariants (for the flopped geometry, or equivalently for anti-M2-branes). The BPS generating functions read 𝒡conifold𝑛=π‘›βˆ’1𝑗=1ξƒ©π‘ž1βˆ’π‘—π‘ π‘„ξƒͺ𝑗.(2.26) The corresponding statistical models were shown in [16, 18, 21] to correspond to finite pyramids with π‘›βˆ’1 stones in the top row, as shown in Figure 3. In this case, the generating functions of such partitions are equal to 𝑍pyramidπ‘›ξ€·π‘ž0,π‘ž1ξ€Έ=π‘›βˆ’1𝑗=1ξ‚€1+π‘ž0π‘›βˆ’π‘—π‘ž1π‘›βˆ’π‘—βˆ’1𝑗.(2.27) The equality 𝒡conifold𝑛≑𝑍pyramid𝑛 arises upon an identification 𝒡conifold𝑛chambersβˆΆπ‘žπ‘ βˆ’1=π‘ž0π‘ž1,𝑄=βˆ’π‘žπ‘›π‘ π‘ž1.(2.28) There are two other sets of chambers characterized by the negative value of the 𝐡-field, for which BPS generating functions are completely analogous to those given above.

Above, we presented just the simplest examples of crystal models. Using fermionic formulation presented below, one can find other crystal models for arbitrary toric geometry without compact four cycles. Let us also mention that those models can be equivalently expressed in terms dimers. In particular, the operation of enlarging the crystal, as in the conifold pyramids, corresponds to so-called dimer shuffling [15]. Dimers are also closely related to a formulation using quivers and associated potentials, which underlies physical derivations in [21, 22].

3. A Little Backgroundβ€”Free Fermions and Matrix Models

In this section, we introduce some mathematical background on which the main results presented in this paper rely. In Section 3.1, we start with a brief presentation of toric Calabi-Yau manifolds and introduce the notation which we use in what follows. In Section 3.2, we introduce free fermion formalism. In Section 3.3, we introduce basics of matrix model formalism. Our presentation is necessarily brief, and for more detailed introduction we recommend many excellent reviews on each of those topics.

3.1. Toric Calabi-Yau Threefolds

Some introductory material on toric Calabi-Yau manifolds, from the perspective relevant for mirror symmetry and topological string thoery, can be found, for example, in [39]. In this section, our presentation is brief and mainly sets up the notation. Toric Calabi-Yau threefolds arise as the quotient of β„‚πœ…+3, possibly with a discrete set of points deleted, by the action of (β„‚βˆ—)πœ… with certain weights. The simplest toric threefold is β„‚3, which corresponds to the trivial choice πœ…=0. The resolved conifold, which we already discussed in Section 2.2, corresponds to πœ…=1 and a choice of weights (1,1,βˆ’1,βˆ’1), which represent a local bundle π’ͺ(βˆ’1)βŠ•π’ͺ(βˆ’1)β†’β„™1. The structure of each toric three-fold can be encoded in a two-dimensional diagram built from trivalent vertices. Finite intervals joining two adjacent vertices represent local β„™1 neighborhood inside the manifold. Equivalently, one can consider dual graphs. Examples of toric diagrams and their duals for β„‚3, conifold and resolution of β„‚3/β„€2 singularity are given in Figure 4 (the notation Γ± at each vertex will be explained in what follows).

A closed loop in a toric diagram represents a compact four cycles in the geometry. As follows from the reasoning in Section 2.1, in the context of BPS counting, we are forced to restrict considerations to manifolds which do not have such four cycles. Apart from a few special cases, there is an infinite class of such geometries whose dual diagrams arise from a triangulation, into triangles of area 1/2, of a long rectangle or a strip of height 1. A toric diagram arises as a dual graph to such a triangulation. From each vertex in such a toric diagram, one vertical line extends to infinity and crosses either the upper or the lower edge of the strip. Two such consecutive lines can emanate either in the same or in the opposite direction, respectively, when they are the endpoints of an interval representing β„™1 with local π’ͺ(βˆ’2)βŠ•π’ͺ or π’ͺ(βˆ’1)βŠ•π’ͺ(βˆ’1) neighborhood. An example of a generic diagram of this kind is shown in Figure 8.

Let us denote independent β„™1's, starting from the left end of the strip, from 1 to 𝑁, and introduce corresponding KΓ€hler parameters 𝑄𝑖=π‘’βˆ’π‘‡π‘–, 𝑖=1,…,𝑁. Moreover, to each toric vertex we associate a type 𝑑𝑖=Β±1, so that 𝑑𝑖+1=𝑑𝑖 if the local neighborhood of β„™1 (represented by an interval between vertices 𝑖 and 𝑖+1) is π’ͺ(βˆ’2)βŠ•π’ͺ; if this neighborhood is of π’ͺ(βˆ’1)βŠ•π’ͺ(βˆ’1) type, then 𝑑𝑖+1=βˆ’π‘‘π‘–. The type of the first vertex we fix as 𝑑1=+1. In Figures 4 and 8, these types are denoted by βŠ• and βŠ–. The types 𝑑𝑖 will be used much in the construction of fermionic states in Section 4.2.

As explained in Section 2.1, the BPS generating functions can be expressed in terms (the instanton part) of topological string amplitudes. For the above class of geometries, arising from a triangulation of a strip, these amplitudes read 𝒡topξ€·π‘„π‘–ξ€Έξ€·π‘ž=π‘€π‘ ξ€Έβˆž(𝑁+1)/2𝑙=11≀𝑖<𝑗≀𝑁+1ξ€·1βˆ’π‘žπ‘™π‘ ξ€·π‘„π‘–π‘„π‘–+1β‹―π‘„π‘—βˆ’1ξ€Έξ€Έβˆ’(𝑑𝑖𝑑𝑗)𝑙.(3.1)

3.2. Free Fermion Formalism

Formalism of free fermions in two dimensions is well known [43, 44] and ubiquitous in literature on topological strings and crystal melting [7, 15, 15, 37, 45]. The main purpose of this section is therefore to set up the notation which we will follow in the remaining parts of this paper.

The states in the free fermion Fock space are created by the (anticommuting) modes of the fermion field ξ“πœ“(𝑧)=π‘˜βˆˆβ„€πœ“π‘˜+1/2π‘§βˆ’π‘˜βˆ’1,πœ“βˆ—(𝑧)=π‘˜βˆˆβ„€πœ“βˆ—π‘˜+1/2π‘§βˆ’π‘˜βˆ’1,ξ€½πœ“π‘˜+(1/2),πœ“βˆ—βˆ’π‘™βˆ’1/2ξ€Ύ=π›Ώπ‘˜,𝑙(3.2) on the vacuum state |0⟩. There is one-to-one map between such fermionic states ||πœ‡βŸ©=𝑑𝑖=1πœ“βˆ—βˆ’π‘Žπ‘–βˆ’1/2πœ“βˆ’π‘π‘–βˆ’1/2||0⟩,withπ‘Žπ‘–=πœ‡π‘–βˆ’π‘–,𝑏𝑖=πœ‡π‘‘π‘–βˆ’π‘–,(3.3) and two-dimensional partitions πœ‡=(πœ‡1,πœ‡2,…,πœ‡π‘™), as shown in Figure 5. The modes π›Όπ‘š of the bosonized field πœ•πœ™=βˆΆπœ“(𝑧)πœ“βˆ—(𝑧)∢ satisfy the Heisenberg algebra [π›Όπ‘š,π›Όβˆ’π‘›]=π‘›π›Ώπ‘š,𝑛.

We introduce vertex operators Γ±(π‘₯)=π‘’βˆ‘π‘›>0(π‘₯𝑛/𝑛)𝛼±𝑛,Ξ“ξ…žΒ±(π‘₯)=π‘’βˆ‘π‘›>0((βˆ’1)π‘›βˆ’1π‘₯𝑛/𝑛)𝛼±𝑛,(3.4) which act on fermionic states |πœ‡βŸ© corresponding to partitions πœ‡ as [15, 43, 44] Ξ“βˆ’(||π‘₯)πœ‡βŸ©=πœ†β‰»πœ‡π‘₯|πœ†|βˆ’|πœ‡|||πœ†βŸ©,Ξ“+(||π‘₯)πœ‡βŸ©=πœ†β‰Ίπœ‡π‘₯|πœ‡|βˆ’|πœ†|||Ξ“πœ†βŸ©,(3.5)ξ…žβˆ’||(π‘₯)πœ‡βŸ©=πœ†π‘‘β‰»πœ‡π‘‘π‘₯|πœ†|βˆ’|πœ‡|||πœ†βŸ©,Ξ“ξ…ž+||(π‘₯)πœ‡βŸ©=πœ†π‘‘β‰Ίπœ‡π‘‘π‘₯|πœ‡|βˆ’|πœ†|||πœ†βŸ©.(3.6) The interlacing relation β‰Ί between partitions is defined as πœ†β‰»πœ‡βŸΊπœ†1β‰₯πœ‡1β‰₯πœ†2β‰₯πœ‡2β‰₯πœ†3β‰₯….(3.7)

The operator Ξ“β€² is the inverse of Ξ“ with negative argument. These operators satisfy commutation relations Ξ“+(π‘₯)Ξ“βˆ’1(𝑦)=Ξ“1βˆ’π‘₯π‘¦βˆ’(𝑦)Ξ“+Ξ“(π‘₯),(3.8)ξ…ž+(π‘₯)Ξ“ξ…žβˆ’1(𝑦)=Ξ“1βˆ’π‘₯π‘¦ξ…žβˆ’(𝑦)Ξ“ξ…ž+Ξ“(π‘₯),(3.9)ξ…ž+(π‘₯)Ξ“βˆ’(𝑦)=(1+π‘₯𝑦)Ξ“βˆ’(𝑦)Ξ“ξ…ž+Ξ“(π‘₯),(3.10)+(π‘₯)Ξ“ξ…žβˆ’(𝑦)=(1+π‘₯𝑦)Ξ“ξ…žβˆ’(𝑦)Ξ“+(π‘₯).(3.11)

We also introduce various colors π‘žπ‘” and the corresponding operators 𝑄𝑔 (a hat is to distinguish them from KΓ€hler parameters 𝑄𝑖) 𝑄𝑔||πœ†βŸ©=π‘žπ‘”|πœ†|||πœ†βŸ©.(3.12) These operators commute with vertex operators up to rescaling of their arguments Ξ“+𝑄(π‘₯)𝑔=𝑄𝑔Γ+ξ€·π‘₯π‘žπ‘”ξ€Έ,Ξ“ξ…ž+𝑄(π‘₯)𝑔=ξπ‘„π‘”Ξ“ξ…ž+ξ€·π‘₯π‘žπ‘”ξ€Έ,𝑄(3.13)π‘”Ξ“βˆ’(π‘₯)=Ξ“βˆ’ξ€·π‘₯π‘žπ‘”ξ€Έξπ‘„π‘”,ξπ‘„π‘”Ξ“ξ…žβˆ’(π‘₯)=Ξ“ξ…žβˆ’ξ€·π‘₯π‘žπ‘”ξ€Έξπ‘„π‘”.(3.14)

3.3. Matrix Models

In matrix model theory, or theory of random matrices, one is interested in properties of various ensembles of matrices. Excellent reviews of random matrix theory can be found for example in [46] or, in particular in the context of topological string theory, in [47]. In matrix model theory, one typically considers partition functions of the form ξ€œξ‘π‘=π’Ÿπ‘ˆπ›Όπ‘’βˆ’(1/𝑔𝑠)Tr𝑉(π‘ˆ),(3.15) where 𝑉=𝑉(π‘ˆ) is a matrix potential and π’Ÿπ‘ˆ is a measure over a set of matrices of interest π‘ˆ of size 𝑁. Typically it is not possible to perform the above integral; however, special techniques allow to determine its formal 1/𝑁 expansion. These techniques culminated with the formalism of the topological expansion of Eynard and Orantin [48] which, in principle, allows to determine entire 1/𝑁 expansion of the partition function recursively. This solution is determined by the behavior of matrix eigenvalues, whose distribution among the minima of the potential, in the continuum limit, determines one-dimensional complex curve, so-called spectral curve. The spectral curve is also encoded in the leading 1/𝑁 expansion of the so-called resolvent, which is defined as the expectation value πœ”(π‘₯)=⟨Tr(1/(π‘₯βˆ’π‘ˆ))⟩ computed with respect to the measure (3.15).

In the context of BPS counting and topological strings, unitary ensembles of matrices of infinite size arise. In this case, the matrix model simplifies to the integral over eigenvalues 𝑒𝛼, with a measure which takes form of the unitary Vandermonde determinant ξ‘π’Ÿπ‘ˆ=𝛼𝑑𝑒𝛼𝛼<𝛽||π‘§π›Όβˆ’π‘§π›½||2,𝑧𝛼=𝑒𝑖𝑒𝛼.(3.16) The issue of infinite matrices is a little subtle; however, it can be taken care of by considering matrices of large but finite size 𝑁, and subsequently taking π‘β†’βˆž limit. For finite 𝑁, one can find the resolvent, and in consequence the spectral curve, using a standard technique of so-called Migdal integral. This requires redefining 𝑉 to the standard Vandermonde form [29, 47], as well as introducing 't Hooft coupling π‘‡π‘‰βŸΆπ‘‰+𝑇log𝑧,𝑇=𝑁𝑔𝑠.(3.17) The form of the Migdal integral depends on the number of cuts into which eigenvalues condense in large 𝑁 limit, and this number of cuts determines the genus of the spectral curve. In our context, only single-cut situations will arise, for which the spectral curve has genus zero. In this case, the Migdal integral determines the resolvent as 1πœ”(𝑝)=ξ€Ÿ2π‘‡π‘‘π‘§πœ•2πœ‹π‘–π‘§π‘‰(𝑧)βˆšπ‘βˆ’π‘§(π‘βˆ’π‘Ž)(π‘βˆ’π‘)√(π‘§βˆ’π‘Ž)(π‘§βˆ’π‘),(3.18) so that the integration contour encircles counter-clockwise the endpoints of the cut π‘Ž and 𝑏. A proper asymptotic behavior of the resolvent is imposed by the condition limπ‘β†’βˆž1πœ”(𝑝)=𝑝.(3.19) Then the spectral curve is determined as a surface on which the resolvent is unambiguously defined, that is, it is given by an (exponential) rational equation automatically satisfied by 𝑝 and πœ”(𝑝). There is also an important consistency condition for the resolvent: when computed on the opposite sides of the cut πœ”(𝑝)Β±, it is related to the potential as πœ”+(𝑝)+πœ”βˆ’πœ•(𝑝)=𝑝𝑉(𝑝)𝑇.(3.20) On the other hand, a difference of these values of the resolvent on both sides of the cut provides eigenvalue density 𝜌(𝑝)=πœ”+(𝑝)βˆ’πœ”βˆ’(𝑝).(3.21)

It has been observed in several contexts that topological strings on toric manifolds can be related to matrix models, whose spectral curves take form of the so-called mirror curves. Mirror curves arise for manifolds which are mirror to toric Calabi-Yau manifolds [39, 45]. For toric manifolds, their mirror manifolds are determined by the following equation embedded in four-dimensional complex space: 𝑧1𝑧2=𝐻(π‘₯,𝑦).(3.22) The mirror curve is the zero locus of 𝐻(π‘₯,𝑦), that is, it is given as 𝐻(π‘₯,𝑦)=0. More precisely, π‘₯,𝑦 are β„‚βˆ— variables, and it is often convenient to represent them in the exponential form π‘₯=𝑒𝑒, 𝑦=𝑒𝑣, with 𝑒,π‘£βˆˆβ„‚. For example, for β„‚3 and the conifold they take the following form: 𝐻ℂ3(π‘₯,𝑦)=π‘₯+𝑦+π‘₯𝑦=0,𝐻conifold(π‘₯,𝑦)=π‘₯+𝑦+π‘₯𝑦+𝑄π‘₯2=0,(3.23) where 𝑄 encodes the KΓ€hler parameter of the conifold. Schematically mirror curves arise from thickening edges of the toric graphs, as shown in Figure 6.

One of the first relations between topological strings for toric manifolds and matrix models was encountered in [49, 50], where it was shown that the spectral curve of a unitary matrix model with a Gaussian (i.e., quadratic) potential agrees with the above mirror curve 𝐻conifold(π‘₯,𝑦)=0 in (3.23), with ’t Hooft coupling 𝑇=𝑔𝑠𝑁 encoded in 𝑄=π‘’βˆ’π‘‡. At the same time, it was shown that the matrix model partition function reproduces the topological string partition function. More recently these ideas became important in view of the remodeling conjecture [51], which states that the solution to loop equations in the form found by Eynard and Orantin [48], applied to the mirror curve, reproduces topological string partition functions. The method of [48] works for arbitrary curves, not necessarily originating from matrix models. Nonetheless, it is indeed possible to construct matrix models whose partition functions do reproduce topological string amplitudes, and whose spectral curves coincide with appropriate mirror curves [29, 30, 52–56].

One of our aims is to provide matrix model interpretation of BPS counting. It is natural to expect such an interpretation in view of an intimate relation between BPS counting and topological string theory discussed in Section 2.1, and the above-mentioned relations between topological strings and matrix models. As we will see in what follows, there are indeed unitary matrix models which naturally arise in the context of BPS counting and its fermionic formulation. Among the others, our task will be to analyze them using the above-mentioned Migdal method.

4. Fermionic Formulation of BPS Counting Functions

Having introduced all the ingredients above, we are now ready to present fermionic formulation of BPS counting. To start with, in Section 4.1 we present the idea of such a formulation in the simplest example of β„‚3. In Section 4.2, we introduce a general fermionic formalism, and in Section 4.3 we provide its crystal interpretation. We illustrate the use of our formalism in Section 4.4 revisiting β„‚3 example, as well as in explicit case of β„‚3/℀𝑁, and conifold geometry.

4.1. The Idea and β„‚3 Example

As explained in Section 2.2, the generating function of bound states of D0-branes to a single D6-brane is given by the MacMahon function, and the corresponding crystal model takes form of the counting of plane partitions [7]. Let us slice each such plane partition by a set of parallel planes, as shown in Figure 7. In this way on each slice, we obtain a two-dimensional partition πœ‡, and it is not hard to see that each two neighboring partitions satisfy the interlacing condition (3.7). Recalling that such a condition arises if we apply Γ±(1) operators (3.5) to partition states, we conclude that a set of all plane partitions can be built, slice by slice, by acting with infinite sequence of Γ±(1) on the vacuum. To count each slice πœ‡ with appropriate weight π‘ž|πœ‡| we also need to apply weight operator 𝑄 defined in (3.12). Therefore, the generating function of plane partitions can be represented as followsΩ𝑍=+βˆ£Ξ©βˆ’ξ¬||β‹―ξβ‰‘βŸ¨0𝑄Γ+(1)𝑄Γ+(1)𝑄Γ+(1)βˆ£π‘„Ξ“βˆ’ξ(1)π‘„Ξ“βˆ’ξ(1)π‘„Ξ“βˆ’ξ||||(1)𝑄⋯0⟩=⟨0β‹―Ξ“+ξ€·π‘ž2ξ€ΈΞ“+(π‘ž)Ξ“+(1)Ξ“βˆ’(π‘ž)Ξ“βˆ’ξ€·π‘ž2ξ€ΈΞ“βˆ’ξ€·π‘ž3ξ€Έβ‹―||=0βŸ©βˆžξ‘π‘™1,𝑙2=111βˆ’π‘žπ‘™1+𝑙2βˆ’1=𝑀(π‘ž).(4.1) In the first line, we implicitly introduced two states ⟨Ω+| and |Ξ©βˆ’βŸ©, defined by an infinite sequence of Ξ“+ (resp. Ξ“βˆ’) operators, interlaced with weight operators 𝑄 and acting on the vacuum. To confirm that this correlator indeed reproduces the MacMahon function, the second line can be reduced to the final infinite product using commutation relations (3.8) and (3.14). We can also represent insertions of Γ±(1) operators graphically by arrows, so that the above computation can be represented as in Figure 7(b).

In what follows, we present a formalism which allows to generalize this computation to a large class of chambers, for arbitrary toric geometry without compact four cycles.

4.2. Toric Geometry and Quantization

We wish to reformulate BPS counting in the fermionic language in a way in which we associate to each toric manifold a fermionic state, such that the BPS generating function can be expressed as an overlap of two such states, generalizing β„‚3 case (4.1). At the same time, the construction of such a fermionic state is supposed to encode the structure of the underlying crystal model (generalizing plane partitions in Figure 7). An important difference between β„‚3 and other geometries is the existence of many KΓ€hler moduli and correspondingly many chambers, for which BPS generating functions change according to wall-crossing formulas. To take care of these changes in the fermionic formalism, we need to introduce special wall-crossing operators.

4.2.1. Toric Geometry and Fermionic Operators

In what follows we use the notation introduced in Section 3.1; in particular to each vertex of the toric diagram we associate its type 𝑑𝑖=Β±1, see also Figure 8. We start with a construction of fermionic states associated to a given toric Calabi-Yau manifold (without compact four-cycles). First we need to introduce several operators which are building blocks of such states. The structure of these operators is encoded in the toric diagram of a given manifold. Namely, these operators are given by a string of 𝑁+1 vertex operators Γ𝑑𝑖±(π‘₯) (defined in (3.4)) which are associated to the vertices of the toric diagram; the type 𝑑𝑖 determines the type of a vertex operator as Γ𝑑𝑖±=+1(π‘₯)=Γ±(π‘₯),Γ𝑑𝑖±=βˆ’1(π‘₯)=Ξ“ξ…žΒ±(π‘₯).(4.2) In addition the string of operators Γ𝑑𝑖±(π‘₯) is interlaced with 𝑁+1 operators 𝑄𝑖 representing colors π‘žπ‘–, for 𝑖=0,1,…,𝑁. Operators 𝑄1𝑄,…,𝑁 are associated to β„™1 in the toric diagram, and there is an additional 𝑄0. We also define 𝑄𝑄=0𝑄1⋯𝑄𝑁,π‘ž=π‘ž0π‘ž1β‹―π‘žπ‘.(4.3) Therefore, the upper indices of Γ𝑑𝑖±(π‘₯) and a choice of colors of the operators which we introduce below are specified by the data of a given toric manifold. As we will see, a sequence of lower indices Β± is determined by the chamber we are going to consider.

Now we can associate several operators to a given toric manifold. Firstly, we define 𝐴±(π‘₯)=Γ𝑑1±𝑄(π‘₯)1Γ𝑑2±𝑄(π‘₯)2⋯Γ𝑑𝑁±𝑄(π‘₯)𝑁Γ𝑑𝑁+1±𝑄(π‘₯)0.(4.4) Commuting all 𝑄𝑖's using (3.14), we also define the following operators: 𝐴+𝑄(π‘₯)=βˆ’1𝐴+(π‘₯)=Γ𝑑1+(π‘₯π‘ž)Γ𝑑2+ξ‚΅π‘₯π‘žπ‘ž1Γ𝑑3+ξ‚΅π‘₯π‘žπ‘ž1π‘ž2⋯Γ𝑑𝑁+1+ξ‚΅π‘₯π‘žπ‘ž1π‘ž2β‹―π‘žπ‘ξ‚Ά,π΄βˆ’(π‘₯)=π΄βˆ’(𝑄π‘₯)βˆ’1=Γ𝑑1βˆ’(π‘₯)Γ𝑑2βˆ’ξ€·π‘₯π‘ž1Γ𝑑3βˆ’ξ€·π‘₯π‘ž1π‘ž2⋯Γ𝑑𝑁+1βˆ’ξ€·π‘₯π‘ž1π‘ž2π‘žπ‘ξ€Έ.(4.5)

In addition, we define the above-mentioned wall-crossing operators π‘Šπ‘ξ‚€Ξ“(π‘₯)=𝑑1βˆ’ξπ‘„(π‘₯)1Γ𝑑2βˆ’ξπ‘„(π‘₯)2β‹―Ξ“π‘‘π‘βˆ’ξπ‘„(π‘₯)𝑝Γ𝑑𝑝+1+𝑄(π‘₯)𝑝+1⋯Γ𝑑𝑁+𝑄(π‘₯)𝑁Γ𝑑𝑁+1+𝑄(π‘₯)0,π‘Šξ…žπ‘ξ‚€Ξ“(π‘₯)=𝑑1+𝑄(π‘₯)1Γ𝑑2+𝑄(π‘₯)2⋯Γ𝑑𝑝+𝑄(π‘₯)𝑝Γ𝑑𝑝+1βˆ’ξπ‘„(π‘₯)𝑝+1β‹―Ξ“π‘‘π‘βˆ’ξπ‘„(π‘₯)𝑁Γ𝑑𝑁+1βˆ’ξπ‘„(π‘₯)0.(4.6) Here the order of Ξ“ and Ξ“ξ…ž is the same as for 𝐴± operators, and the difference is that now there are subscripts βˆ“ on first 𝑝 operators and Β± on the remaining ones.

We often use a simplified notation when the argument of the above operators is π‘₯=1𝐴±≑𝐴±(1),𝐴±≑𝐴±(1),π‘Šπ‘β‰‘π‘Šπ‘(1),π‘Šξ…žπ‘β‰‘π‘Šξ…žπ‘(1).(4.7)

4.2.2. Fermionic Formulation and Quantization

Above we associated operators 𝐴± to each toric geometry with a strip-like toric diagram. From these operators, we can build the following states in the Hilbert space of a free fermion β„‹: ||Ξ©Β±ξ¬βˆˆβ„‹,(4.8) which we define as follows: Ω+||||β‹―=⟨0𝐴+(1)𝐴+(1)𝐴+||(1)=⟨0⋯𝐴+ξ€·π‘ž2𝐴+(π‘ž)𝐴+||Ξ©(1),βˆ’βŸ©=π΄βˆ’(1)π΄βˆ’(1)π΄βˆ’(||1)β‹―0⟩=π΄βˆ’(1)π΄βˆ’(π‘ž)π΄βˆ’ξ€·π‘ž2ξ€Έβ‹―||0⟩.(4.9) These states encode the full instanton part of the topological string amplitudes. Namely, as shown in [16], Ω𝑍=+βˆ£Ξ©βˆ’ξ¬(4.10) is equal to the BPS partition function 𝒡 in the noncommutative chamber ||𝒡𝑍=𝒡≑top||2≑𝒡top𝑄𝑖𝒡topξ€·π‘„π‘–βˆ’1ξ€Έ,(4.11) where 𝒡top(𝑄𝑖) is given in (3.1). The above equality holds under the following identification between π‘žπ‘– parameters (which enter the definition of |Ω±⟩) and physical parameters 𝑄𝑖=π‘’βˆ’π‘‡π‘– and π‘žπ‘ =π‘’βˆ’π‘”π‘ : π‘žπ‘–=𝑑𝑖𝑑𝑖+1𝑄𝑖,π‘žπ‘ =π‘žβ‰‘π‘ž0π‘ž1β‹―π‘žπ‘.(4.12) We will provide a proof of (4.10) in Section 6.1.1 in a more general setting of refined invariants.

The states |Ω±⟩ have nontrivial structure and encode the information about the noncommutative chamber. It turns out that the fermionic vacuum |0⟩ itself also encodes some interesting information. We recall that there is another extreme chamber representing just a single BPS state represented by the D6-brane with no other branes bound to it. This multiplicity 1 can be understood as 𝒡=𝑍=⟨0∣0⟩=1,(4.13) and as we will see below, starting from this expression we can use wall-crossing operators to construct BPS generating functions in an infinite family of other chambers.

4.2.3. Other Chambers and Wall-Crossing Operators

In the previous, section we associated to toric manifolds the states |Ω±⟩, whose overlap reproduces the BPS generating function in the noncommutative chamber (4.10). Now we wish to extend this formalism to other chambers. As discussed in Section 2.1, in a given chamber, the allowed bound states we wish to count must have positive central charge (2.1) 1𝑍(𝑅,𝐡)=𝑅(𝑛+𝛽⋅𝐡)>0.(4.14) Firstly, the information about 𝑅 and 𝐡 must be encoded in the fermionic states which we wish to construct. It turns out that the choice of positive or negative 𝑅 is encoded in the choice of the ground state ||Ω𝑅>0⟢±||,𝑅<0⟢0⟩,(4.15) which generalizes the extreme cases (4.10) and (4.13).

On the other hand, the value of the field 𝐡 is encoded in the insertion of additional wall-crossing operators, such as those defined in (4.6). In particular, these two types of operators are sufficient if we wish to consider only these chambers, which correspond to a flux of the 𝐡-field through only one, but arbitrary β„™1 in the manifold. For simplicity below, we consider only this set of chambers. Denoting this β„™1 as 𝑝, it can be shown that insertion of 𝑛 copies of operators π‘Šπ‘ or π‘Šξ…žπ‘ creates, respectively, 𝑛 positive or negative quanta of the flux through 𝑝'th β„™1.

Therefore, schematically, the generating functions in chambers with 𝑅>0 read 𝑍𝑛=Ω+||ξ‚€π‘Šξ‚π‘›||Ξ©βˆ’βŸ©,(4.16) and those with 𝑅<0 read 𝑍𝑛||ξ‚€=⟨0π‘Šξ‚π‘›||0⟩,(4.17) with appropriate form of wall-crossing operators. More precisely, depending on the signs of 𝑅 and 𝐡, we need to consider four possible situations, which we present below. The proofs of all statements below, corresponding to these four situations, can be found in [16].

(i) Chambers with 𝑅<0, 𝐡>0
Consider a chamber characterized by positive 𝑅 and positive 𝐡-field through 𝑝'th two-cycle ][𝑅<0,π΅βˆˆπ‘›βˆ’1,𝑛,for1β‰€π‘›βˆˆβ„€.(4.18) The BPS partition function in this chamber contains only those factors which include 𝑄𝑝 and it reads 𝒡𝑛|𝑝=π‘›βˆ’1𝑝𝑖=1𝑠=1𝑁+1ξ‘π‘Ÿ=𝑝+1ξ‚΅π‘ž1βˆ’π‘–π‘ π‘„π‘ π‘„π‘ +1β‹―π‘„π‘Ÿβˆ’1ξ‚Άβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘–.(4.19) This can be expressed as the expectation value of 𝑛 wall-crossing operators π‘Šπ‘ξ‚π‘π‘›|𝑝||ξ‚€=⟨0π‘Šπ‘ξ‚π‘›||𝒡0⟩=𝑛|𝑝,(4.20) under the following identification of variables: 𝑄𝑝=𝑑𝑝𝑑𝑝+1ξ€Έπ‘žπ‘π‘žπ‘›π‘ ,𝑄𝑖=𝑑𝑖𝑑𝑖+1ξ€Έπ‘žπ‘–for𝑖≠𝑝,π‘žπ‘ =1π‘ž.(4.21) A special case of this result is the trivial generating function (4.13) representing a single D6-brane.

(ii) Chambers with 𝑅>0, 𝐡>0
In the second case, we consider the positive value of 𝑅 and the positive flux through 𝑝'th β„™1][𝑅>0,π΅βˆˆπ‘›,𝑛+1,for0β‰€π‘›βˆˆβ„€.(4.22) Denote the BPS partition function in this chamber by 𝒡𝑛|𝑝. We find that the expectation value of 𝑛 wall-crossing operators π‘Šπ‘ in the background of |Ω⟩ has the form 𝑍𝑛|𝑝=Ω+||ξ‚€π‘Šπ‘ξ‚π‘›||Ξ©βˆ’βŸ©=𝑀(1,π‘ž)𝑁+1𝑍(0)𝑛|𝑝𝑍(1)𝑛|𝑝𝑍(2)𝑛|𝑝,(4.23) where 𝑍(0)𝑛|𝑝 does not contain any factors (π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1)Β±1 which would include π‘žπ‘, while 𝑍(1)𝑛|𝑝 contains all factors π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1 which do include π‘žπ‘, and 𝑍(2)𝑛|𝑝 contains all factors (π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1)βˆ’1 which also include π‘žπ‘: 𝑍(0)𝑛|𝑝=βˆžξ‘π‘™=1ξ‘π‘βˆ‰π‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ‚Άβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™ξ€·ξ€·π‘‘1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ€Έβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™,𝑍(1)𝑛|𝑝=βˆžξ‘π‘™=1ξ‘π‘βˆˆπ‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™+π‘›π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ€Έβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™,𝑍(2)𝑛|𝑝=βˆžξ‘π‘™=𝑛+1ξ‘π‘βˆˆπ‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™βˆ’π‘›π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ‚Άβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™.(4.24) We see that the identification of variables 𝑄𝑝=𝑑𝑝𝑑𝑝+1ξ€Έπ‘žπ‘π‘žπ‘›π‘ ,𝑄𝑖=𝑑𝑖𝑑𝑖+1ξ€Έπ‘žπ‘–for𝑖≠𝑝,π‘žπ‘ =π‘ž(4.25) reproduces the BPS partition function 𝒡𝑛|𝑝=𝑍𝑛|𝑝.(4.26) When no wall-crossing operator is inserted the change of variables reduces to (4.12) and we get the noncommutative Donaldson-Thomas partition function (4.11), 𝑍0|𝑝=𝒡.

(iii) Chambers with R<0, B<0
Now we consider negative 𝑅 and negative 𝐡-field ][𝑅<0,π΅βˆˆβˆ’π‘›βˆ’1,βˆ’π‘›for0β‰€π‘›βˆˆβ„€.(4.27) For such a chamber the BPS partition function reads ξ‚π’΅ξ…žπ‘›|𝑝=𝑛𝑝𝑖=1𝑠=1𝑁+1ξ‘π‘Ÿ=𝑝+1ξ€·1βˆ’π‘žπ‘–π‘ π‘„π‘ π‘„π‘ +1β‹―π‘„π‘Ÿβˆ’1ξ€Έβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘–.(4.28) Now we find the expectation value of 𝑛 wall-crossing operators π‘Šξ…žπ‘ is equal to ξ‚π‘ξ…žπ‘›|𝑝||ξ‚€=⟨0π‘Šξ…žπ‘ξ‚π‘›||𝒡0⟩=ξ…žπ‘›|𝑝,(4.29) under the change of variables 𝑄𝑝=𝑑𝑝𝑑𝑝+1ξ€Έπ‘žπ‘π‘žπ‘ βˆ’π‘›,𝑄𝑖=𝑑𝑖𝑑𝑖+1ξ€Έπ‘žπ‘–for𝑖≠𝑝,π‘žπ‘ =1π‘ž.(4.30) Now an insertion of π‘Šπ‘ has an interpretation of turning on a negative quantum of 𝐡-field, and the redefinition of 𝑄𝑝 can be interpreted as effectively reducing 𝑑𝑝 by one unit of 𝑔𝑠. As already discussed, ξ‚π‘ξ…ž0βˆ£π‘=⟨0∣0⟩=1(4.31) represents a chamber with a single D6-brane and no other branes bound to it.

(iv) Chambers with 𝑅>0, 𝐡<0
In the last case, we consider positive 𝑅 and negative 𝐡][𝑅>0,0>π΅βˆˆβˆ’π‘›,βˆ’π‘›+1,for1β‰€π‘›βˆˆβ„€.(4.32) We denote the BPS partition function in this chamber by π’΅ξ…žπ‘›|𝑝. We find that the expectation value of 𝑛 operators π‘Šξ…žπ‘ in the background of |Ω±⟩ has the form π‘ξ…žπ‘›|𝑝=Ω+||ξ‚€π‘Šξ…žπ‘ξ‚π‘›||Ξ©βˆ’βŸ©=𝑀(1,π‘ž)𝑁+1π‘ξ…ž(0)𝑛|π‘π‘ξ…ž(1)𝑛|π‘π‘ξ…ž(2)𝑛|𝑝,(4.33) where π‘ξ…ž(0)𝑛|𝑝 does not contain any factors (π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1)Β±1 which would include π‘žπ‘, π‘ξ…ž(1)𝑛|𝑝 contains all factors π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1 which do include π‘žπ‘, and π‘ξ…ž(2)𝑛|𝑝 contains all factors (π‘žπ‘ β‹―π‘žπ‘Ÿβˆ’1)βˆ’1 which also include π‘žπ‘: π‘ξ…ž(0)𝑛|𝑝=βˆžξ‘π‘™=1ξ‘π‘βˆ‰π‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ‚Άβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™ξ€·ξ€·π‘‘1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ€Έβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™,π‘ξ…ž(1)𝑛|𝑝=βˆžξ‘π‘™=π‘›ξ‘π‘βˆˆπ‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™βˆ’π‘›π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ€Έβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™,π‘ξ…ž(2)𝑛|𝑝=βˆžξ‘π‘™=1ξ‘π‘βˆˆπ‘ ,π‘Ÿ+1βŠ‚1,𝑁+1𝑑1βˆ’π‘Ÿπ‘‘π‘ ξ€Έπ‘žπ‘™+π‘›π‘žπ‘ π‘žπ‘ +1β‹―π‘žπ‘Ÿβˆ’1ξ‚Άβˆ’π‘‘π‘Ÿπ‘‘π‘ π‘™.(4.34) Under the change of variables, 𝑄𝑝=𝑑𝑝𝑑𝑝+1ξ€Έπ‘žπ‘π‘žπ‘ βˆ’π‘›βˆ’1,𝑄𝑖=𝑑𝑖𝑑𝑖+1ξ€Έπ‘žπ‘–,for𝑖≠𝑝,π‘žπ‘ =π‘ž.(4.35) This reproduces the BPS partition function π’΅ξ…žπ‘›|𝑝=π‘ξ…žπ‘›|𝑝.(4.36) We note that both π‘ξ…ž1|𝑝 with the above change of variables, as well as 𝑍0|𝑝 given in (4.23) with a different change of variables in (4.25), lead to the same BPS generating function 𝒡 which corresponds to the noncommutative Donaldson-Thomas invariants.

4.3. Crystal Melting Interpretation

In the previous section, we found a free fermion representation of D6-D2-D0 generating functions. The fermionic correlators which reproduce BPS generating functions automatically provide melting crystal interpretation of these functions [16], generalizing models of plane partitions (for β„‚3) or pyramid partitions (for the conifold), presented in Section 2.2. These crystals are also equivalent to those found in [17, 22].

The crystal interpretation is a consequence of the fact that all operators used in the construction of states |Ω±⟩, as well as the wall-crossing operators, are built just from vertex operators Γ± and Γ± with argument 1, and color operators 𝑄𝑖. As follows from (3.5) and (3.6), insertion of these vertex operators is equivalent to the insertion of two-dimensional partitions satisfying interlacing, or transposed interlacing conditions. An infinite sequence of such interlacing partitions effectively builds up a three-dimensional crystal. A relative position of two adjacent slices is determined by a type of two corresponding vertex operators. On the other hand, insertions of color operators have an interpretation of coloring the crystal. The colors 𝑄𝑖 appear in the same order in each composite operator, so these colors are always repeated periodically in the full correlators. Therefore, three-dimensional crystals are built of interlacing, periodically colored slices.

To get more insight about a geometric structure of a crystal, it is convenient to introduce the following graphical representation. We associate various arrows to the vertex operators, as shown in Figure 9. These arrows follow the order of the vertex operators in the fermionic correlators and are drawn from left to right, or up to down (either of these directions is independent of the orientation of the arrow). Following the order of the vertex operators in a given correlator, and drawing a new arrow at the end of the previous one, produces a zig-zag path which represents a shape of the crystal. The coloring of the crystal is taken care of by keeping track of the order of 𝑄𝑖 operators, and by drawing at the endpoint of each arrow a (dashed) line, rotated by 45Β°, colored according to 𝑄𝑖 which we come across. These lines represent two-dimensional slices in appropriate colors. In this way, the corners of two-dimensional partitions arising from slicing of the crystal are located at the end-points of the arrows. The orientation of arrows represents the interlacing condition (i.e., arrows point from a larger to smaller partition). The interlacing pattern between two consecutive slices corresponds to the types of two consecutive arrows. Finally, the points from which two arrows point outwards represent those stones in the crystal, which can be removed from the initial, full crystal configuration. In fermionic correlators, these points correspond to Γ𝑑𝑖+ followed by Ξ“π‘‘π‘—βˆ’ operators. We illustrate this graphical construction in a few examples in the next section.

4.4. Examples
4.4.1. Revisiting β„‚3

Let us reconsider β„‚3 geometry which motivated our discussion in Section 4.1. In this case, the dual toric diagram consists just of one triangle, see Figure 10(a), so there is just one vertex and only one color 𝑄0≑𝑄, and the operators (4.4) take form 𝐴±=Γ±(1)𝑄.(4.37) In consequence, the BPS partition function (4.10) takes exactly the form (4.1).

The crystal structure can be read off from a sequence of arrows associated to 𝐴± operators, following the rules in Figure 9. This gives rise to the crystal shown in Figure 10(b). This is the same crystal as in Figure 7, which represents plane partitions, however, now seen from the opposite side.

4.4.2. Orbifolds β„‚3/℀𝑁+1

Now we consider the resolution of β„‚3/℀𝑁+1 orbifold. In this case, the toric diagram takes form of a triangle of area (𝑁+1)/2, see Figure 11(a). There are 𝑁 independent β„™1's and 𝑁+1 vertices of the same 𝑑𝑖=+1, and operators in (4.4) take the form 𝐴±=Γ±𝑄(1)1Γ±𝑄(1)2⋯Γ±𝑄(1)𝑁Γ±𝑄(1)0.(4.38) In the noncommutative chamber, the corresponding crystal consists of plane partitions, however, with slices colored periodically in 𝑁+1 colors. The partition function in the noncommutative chamber is given by (4.10).

If we turn on an arbitrary 𝐡-field through a fixed β„™1, the structure of wall-crossing operators gives rise to modified containers, see for example Figure 11(b). In particular, enlarging the 𝐡-field by one unit adds one more yellow corner to the crystal.

The crystals corresponding to 𝑅<0 are also easy to find. In the extreme chamber, we get a trivial (empty) crystal, representing a single D6-brane (4.13). Adding wall-crossing operators results in a crystal with several corners, finite along two axis (and extending infinitely along the third axis), as shown in Figure 11(c).

4.4.3. Resolved Conifold

We already presented pyramid crystals for the conifold in Section 2.2. They arise from our formalism as follows. The dual toric diagram for the conifold, see Figure 12(a), consists of two triangles and encodes a single (𝑁=1) β„™1. Two vertices of the toric diagram correspond to two colors 𝑄1 and 𝑄0, so that 𝑄𝑄=1𝑄0,π‘ž=π‘ž1π‘ž0.(4.39) The operators (4.4) in this case read 𝐴±(π‘₯)=Γ±(π‘₯)𝑄1Ξ“ξ…žΒ±(π‘₯)𝑄0,(4.40) while (4.5) is𝐴+(π‘₯)=Ξ“+(π‘₯π‘ž)Ξ“ξ…ž+ξ‚΅π‘₯π‘žπ‘ž1ξ‚Ά,π΄βˆ’(π‘₯)=Ξ“βˆ’(π‘₯)Ξ“ξ…žβˆ’ξ€·π‘₯π‘ž1ξ€Έ,(4.41) and they satisfy 𝐴+(π‘₯)π΄βˆ’ξ€·(𝑦)=1+π‘₯π‘¦π‘ž/π‘ž1ξ€Έξ€·1+π‘₯π‘¦π‘žπ‘ž1ξ€Έ(1βˆ’π‘₯π‘¦π‘ž)2π΄βˆ’(𝑦)𝐴+(π‘₯).(4.42)

The quantum states (4.9) take form ||Ξ©βˆ’βŸ©=π΄βˆ’(1)π΄βˆ’(π‘ž)π΄βˆ’ξ€·π‘ž2ξ€Έβ‹―||Ω0⟩,+||||=⟨0⋯𝐴+ξ€·π‘ž2𝐴+(π‘ž)𝐴+(1),(4.43) and the wall-insertion operators (4.6) are π‘Š1(π‘₯)=Ξ“βˆ’(π‘₯)𝑄1Ξ“ξ…ž+(π‘₯)𝑄0,π‘Šξ…ž1(π‘₯)=Ξ“+(π‘₯)𝑄1Ξ“ξ…žβˆ’(π‘₯)𝑄0.(4.44)

Therefore, the fermionic correlators take form 𝑍𝑛|1=Ω+||ξ‚€π‘Š1𝑛||Ξ©βˆ’ξ‚π‘βŸ©,𝑛|1||ξ‚€=⟨0π‘Š1𝑛||0⟩.(4.45) and encode generating functions (2.25) and (2.27) introduced in Section 2.2. In the noncommutative chamber, we get the result found first in [14], 𝑍0|1=⟨Ω+βˆ£Ξ©βˆ’βŸ©, while a single D6-brane is encoded in 𝑍=⟨0∣0⟩=1. These crystals are shown in Figures 12(b) and 13.

5. Matrix Models and Open BPS Generating Functions

In this section, we explain how matrix model formalism can be applied to analyze BPS counting functions. In the first part, Section 5.1, we explain how to relate fermionic formalism, derived in the previous section, to matrix model representation. In Section 5.2, we illustrate how to construct matrix models for the closed noncommutative chamber. In Section 5.3, we analyze in detail BPS generating functions for the conifold for all chambers with 𝑅>0, and derive corresponding spectral curves. We discuss how these curves relate to (and generalize) mirror curves, which we find (as we should) in the commutative chamber. In Section 5.4, we reveal that matrix model representation in fact encodes open BPS generating functions, which can be identified with matrix model integrands.

5.1. Matrix Models from Free Fermions

Let us explain how to relate fermionic representation of BPS amplitudes, introduced in Section 4.2, to matrix models. This relies on introducing into fermionic correlators representing BPS generating functions, such as (4.10) or (4.23), a special representation of the identity operator 𝕀. The representation we are interested in also consists of infinite product of vertex operators and arises as follows [29]. Firstly, we can use the representation as a complete set of states 𝕀=|π‘…βŸ©βŸ¨π‘…|, which represent two-dimensional partitions. Using orthogonality relations of π‘ˆ(∞) characters πœ’π‘…, and the fact that these characters are given in terms of Schur functions πœ’π‘…=𝑠𝑅(⃗𝑧) for ⃗𝑧=(𝑧1,𝑧2,𝑧3,…), we can write 𝕀=𝑅||||=ξ“π‘…βŸ©βŸ¨π‘…π‘ƒ,𝑅𝛿𝑃𝑑𝑅𝑑||||=ξ€œξ“π‘ƒβŸ©βŸ¨π‘…π’Ÿπ‘ˆπ‘ƒ,𝑅𝑠𝑃𝑑⃗𝑧𝑠𝑅𝑑||||=ξ€œξƒ©ξ‘βƒ—π‘§π‘ƒβŸ©βŸ¨π‘…π’Ÿπ‘ˆπ›ΌΞ“ξ…žβˆ’ξ€·π‘§π›Όξ€Έ||||0⟩ξƒͺξƒ©βŸ¨0π›ΌΞ“ξ…ž+ξ€·π‘§π›Όβˆ’1ξ€Έξƒͺ.(5.1) When such a representation of the identity operator is introduced into (4.10) or (4.23) (or any other correlator of similar structure), we can commute away Γ𝑑𝑖± operators and get rid of operator expressions. For example, inserting the above identity operator in the string of 𝐴+ operators in (4.16) leads to a matrix model with the unitary measure 𝑍𝑛||=⟨0βˆžξ‘π‘–=π‘˜π΄+||𝕀||(1)π‘˜βˆ’1𝑗=0𝐴+(1)βˆ£π‘Šπ‘›||Ξ©βˆ’βŸ©=ξ€œ||π’Ÿπ‘ˆβŸ¨0βˆžξ‘π‘–=π‘˜+1𝐴+(1)βˆ£π›ΌΞ“ξ…žβˆ’ξ€·π‘§π›Όξ€Έ||||0⟩⟨0π›ΌΞ“ξ…ž+ξ€·π‘§π›Όβˆ’1ξ€Έβˆ£π‘˜ξ‘π‘—=0𝐴+(1)βˆ£π‘Šπ‘›||Ξ©βˆ’βŸ©=π‘“π‘˜π‘›ξ€·π‘ž,π‘„π‘–ξ€Έξ€œξ‘π’Ÿπ‘ˆπ›Όπ‘’βˆ’(1/𝑔𝑠)π‘‰π‘˜π‘›(𝑧𝛼).(5.2) The product over 𝛼 represents distinct eigenvalues 𝑧𝛼. Note that we have inserted 𝕀 at the position π‘˜ in the string of 𝐴+(1) operators. In particular, this affects the form of the resulting potential π‘‰π‘˜π‘›(𝑧). Moreover, apart from matrix integral, we find some overall factors π‘“π‘˜π‘›(π‘ž,𝑄𝑖) which take form of various infinite products. They arise, in a generic chamber, from commutations between Γ± ingredients of wall-crossing operators, and Ξ“βˆ“ ingredients of |Ξ©βˆ“βŸ© states. In the closed noncommutative chamber 𝑛=0, these factors are trivial, 𝑓𝑛=0(π‘ž,𝑄𝑖)=1, and they largely simplify in the commutative chamber π‘›β†’βˆž.

There is a large freedom in choosing the value of k, and it is natural to ask if this choice has some physical interpretation. It was argued in [27] that this is indeed the case, and the choice of π‘˜ is equivalent to the choice of open BPS chamber (open BPS chambers were introduced in Section 2.1). In particular, it turns out that the open generating parameter can be identified with matrix eigenvalues 𝑧𝛼, and the open BPS generating function (2.16) in the open chamber labeled by π‘˜ can be identified with matrix integrand 𝒡openBPS=π‘’βˆ’(1/𝑔𝑠)π‘‰π‘˜π‘›(𝑧).(5.3) Even though the overall factors π‘“π‘˜π‘› in (5.1) may involve closed moduli 𝑄𝑖, they do not involve open moduli 𝑧. In this sense, the matrix integrand is well defined, and up to some simple identification can be identified with open BPS generating function. This identification of parameters amounts to the shift π‘§β†’βˆ’π‘§π‘ž1/2 (to match earlier M-theory convention with half-integer powers of π‘ž, to integer powers of π‘ž in the fermionic formalism), as well as identification of KΓ€hler parameters considered in M-theory derivation with parameters πœ‡π‘– introduced below. We also note that the BPS generating function in (2.16) is determined by the open topological string partition function associated to the external axis of the toric diagram, as in Figure 14. As we will also see, the value of the above integral (5.1) can be related to some more general Calabi-Yau geometry π‘Œ.

5.2. Matrix Models for the Noncommutative Chamber

In this section, we illustrate the relation between BPS counting and matrix models in case of the noncommutative chamber 𝑛=0, and the choice of open chamber also π‘˜=0. This corresponds to the insertion of the identity representation (75) exactly in between |Ω±⟩ states in (4.10). In this 𝑛=0 case no factor π‘“π‘˜π‘› in (5.1) arises, and we obtain matrix models with potentials which can be expressed in terms of the following version of the theta function: Θ(𝑧;π‘ž)=βˆžξ‘π‘—=0ξ€·1+π‘§π‘žπ‘—ξ€Έξ‚΅π‘ž1+𝑗+1𝑧.(5.4) For a general geometry of the form shown in Figure 8, with types of vertices given by 𝑑𝑖, corresponding matrix models take form ξ€œξ‘π‘=π‘‘π‘ˆπ›Όπ‘’βˆ’(1/𝑔𝑠)𝑉(𝑧𝛼)=ξ€œξ‘π‘‘π‘ˆπ›Όπ‘ξ‘π‘™=0Ξ˜ξ€·π‘‘π‘™+1π‘§π›Όξ€·π‘ž1β‹―π‘žπ‘™ξ€Έξ€Έ;π‘žπ‘‘π‘™+1,(5.5) where integral is over unitary matrices of infinite size, 𝑁=∞. Special cases of this result include (i)for β„‚3, the result (5.5) provides a matrix model representation of MacMahon function βˆπ‘=𝑀(1)=βˆžπ‘˜=1(1βˆ’π‘žπ‘˜)βˆ’π‘˜ in terms of a matrix model of the form (5.5) with the integrand π‘’βˆ’(1/𝑔𝑠)𝑉(𝑧)=βˆžξ‘π‘—=0ξ€·1+π‘§π‘žπ‘—ξ€Έξ‚΅π‘ž1+𝑗+1𝑧=Θ(𝑧;π‘ž);(5.6)(ii)for the conifold, we obtain a representation of the pyramid partition generating function (2.23) (with 𝑛=0) in terms of a matrix model with the integrand π‘’βˆ’(1/𝑔𝑠)𝑉(𝑧)=βˆžξ‘π‘—=0ξ€·1+π‘§π‘žπ‘—ξ€Έξ€·1+π‘žπ‘—+1ξ€Έ/𝑧(1+π‘„π‘§π‘žπ‘—)ξ€·ξ€·π‘ž1+𝑗+1=/π‘„π‘§ξ€Έξ€ΈΞ˜(𝑧;π‘ž)Θ(𝑄𝑧;π‘ž);(5.7)(iii)for β„‚3/℀𝑁+1, we have 𝑑𝑝=+1 for all 𝑝 and we find matrix model representation of the BPS generating function in terms of a matrix model with the integrand π‘’βˆ’(1/𝑔𝑠)𝑉(𝑧)=βˆžξ‘π‘—=0ξ€·1+π‘§π‘žπ‘—ξ€Έξ‚΅π‘ž1+𝑗+1π‘§ξ‚Άβ‹―ξ€·ξ€·π‘ž1+1β‹―π‘žπ‘ξ€Έπ‘§π‘žπ‘—ξ€Έξƒ©π‘ž1+𝑗+1ξ€·π‘ž1β‹―π‘žπ‘ξ€Έπ‘§ξƒͺ=𝑁𝑙=0Ξ˜π‘žξ€·ξ€·1β‹―π‘žπ‘™ξ€Έξ€Έ.𝑧;π‘ž(5.8)

5.3. Matrix Model for the Conifold Analysis

In this section, we illustrate how matrix model techniques can be used in the context of models which arise for BPS counting. We focus on the conifold matrix model in arbitrary closed BPS chamber 𝑛, and fixed π‘˜=0. In this case, the result (5.2) takes form (after the redefinition 𝑄=βˆ’π‘ž1π‘žπ‘›) 𝑍𝑛=𝑀(π‘ž)2βˆžξ‘π‘—=1ξ€·1βˆ’π‘„π‘žπ‘—ξ€Έπ‘—ξ€·1βˆ’π‘„βˆ’1π‘žπ‘—+𝑛𝑗+𝑛=π‘“π‘›ξ€œξ‘(π‘ž,𝑄)π‘‘π‘ˆπ›Όβˆžξ‘π‘—=0ξ€·1+π‘§π›Όπ‘žπ‘—+1ξ€Έξ€·1+π‘žπ‘—/𝑧𝛼1+π‘§π›Όπ‘žπ‘—+𝑛+1/𝑄1+π‘žπ‘—π‘„/𝑧𝛼=𝑓𝑛(π‘ž,𝑄)𝑍matrix,(5.9) with π‘“π‘›βˆ(π‘ž,𝑄)=𝑀(π‘ž)βˆžπ‘—=1ξ€·1βˆ’π‘žπ‘›+𝑗/𝑄𝑛𝑀(π‘žπ‘›,π‘ž),(5.10) with MacMahon function 𝑀(π‘ž) defined in (2.17), and with the following generalized MacMahon function 𝑀(𝑧,π‘ž)=βˆžξ‘π‘–=11(1βˆ’π‘§π‘žπ‘–)𝑖.(5.11) In particular, in the noncommutative chamber 𝑓conifold0=1, and in the commutative chamber 𝑓conifoldπ‘›β†’βˆž=𝑀(π‘ž) which represents topological string degree zero contributions. The result (5.9) implies that the value of the matrix model integral (without the prefactor 𝑓𝑛) is equal to 𝑍matrix=𝑍𝑛𝑓𝑛(π‘ž,𝑄)=𝑀(π‘ž)βˆžξ‘π‘—=1ξ€·1βˆ’π‘„π‘žπ‘—ξ€Έπ‘—ξ€·1βˆ’πœ‡π‘žπ‘—ξ€Έπ‘—(1βˆ’πœ‡π‘„π‘žπ‘—)𝑗=ξ€œξ‘π‘‘π‘ˆπ›Όβˆžξ‘π‘—=0ξ€·1+π‘§π›Όπ‘žπ‘—+1ξ€Έξ€·1+π‘žπ‘—/𝑧𝛼1+π‘§π›Όπ‘žπ‘—+𝑛+1/𝑄1+π‘žπ‘—π‘„/𝑧𝛼,(5.12) where πœ‡=π‘žπ‘›/𝑄.

Now we wish to analyze the matrix model 𝑍matrix. We parametrize the 't Hooft coupling and the chamber dependence, respectively, by 𝑇=𝑔𝑠𝑁,𝜏=𝑛𝑔𝑠.(5.13) As our models correspond to π‘ˆ(∞) matrices, ultimately we are interested in the limit π‘‡β†’βˆž,𝑔𝑠=const,𝑄=const,(5.14) for each fixed chamber (i.e., fixed 𝑛 and therefore 𝜏). The noncommutative chamber corresponds to 𝜏=0, while πœβ†’βˆž represents the topological string chamber.

Using the expansion of the quantum dilogarithm logβˆžξ‘π‘–=1ξ€·1βˆ’π‘§π‘žπ‘–ξ€Έ1=βˆ’π‘”π‘ βˆžξ“π‘š=0Li2βˆ’π‘šπ΅(𝑧)π‘šπ‘”π‘šπ‘ π‘š!,(5.15) and the redefinition of the unitary measure (3.17) we find, to the leading order in 𝑔𝑠, the following matrix model potential: π‘‰πœ=𝑇log(𝑧)+Li2(βˆ’π‘§)+Li2ξ‚€βˆ’1π‘§ξ‚βˆ’Li2ξ‚΅βˆ’π‘„π‘§ξ‚Άβˆ’Li2ξ‚΅βˆ’π‘§π‘„π‘’πœξ‚Ά,(5.16) so that πœ•π‘§π‘‰πœ=π‘‡βˆ’log(𝑧+𝑄)+log(1+(𝑧/π‘„π‘’πœ))𝑧.(5.17)

Now we wish to solve the model (5.9) in the small 𝑔𝑠 limit. Firstly, we need to find the resolvent πœ”(𝑝), which can be done using the Migdal integral (3.18), and careful derivation is presented in [29]. As we expect one-cut solution of our model, from the Migdal integral we get an expression in terms of the end-points of this cut π‘Ž and 𝑏. The normalization condition (3.19) imposes two constraints, for terms of order 𝑝0 and π‘βˆ’1 in the resolvent, which take form βˆšβˆšπ‘Ž+π‘„βˆ’π‘+π‘„βˆšπ‘Ž+π‘„π‘’πœβˆ’βˆšπ‘+π‘„π‘’πœ=𝑄1/2𝑒(𝜏+𝑇)/2,√√(π‘Ž+𝑄)π‘βˆ’(𝑏+𝑄)π‘Žβˆš(π‘Ž+π‘„π‘’πœβˆš)π‘βˆ’(𝑏+π‘„π‘’πœ)π‘Ž=𝑄1/2π‘’βˆ’(𝜏+𝑇)/2.(5.18) These constraints can be solved in the exact form, with result π‘Ž=βˆ’1+πœ–2ξ€·(1βˆ’πœ‡)1βˆ’πœ‡πœ–2ξ€Έξ€·+(1βˆ’π‘„)1+πœ‡πœ–2ξ€Έβˆ’2πœ‡ξ€·1βˆ’πœ‡πœ–2ξ€Έ2+2π‘–πœ–ξ€·(1βˆ’π‘„)1βˆ’πœ–2ξ€Έξ€·(1βˆ’πœ‡)1βˆ’π‘„πœ‡πœ–2ξ€Έξ€·1βˆ’πœ‡πœ–2ξ€Έ2,𝑏=βˆ’1+πœ–2ξ€·(1βˆ’πœ‡)1βˆ’πœ‡πœ–2ξ€Έξ€·+(1βˆ’π‘„)1+πœ‡πœ–2ξ€Έβˆ’2πœ‡ξ€·1βˆ’πœ‡πœ–2ξ€Έ2ξ”βˆ’2π‘–πœ–ξ€·(1βˆ’π‘„)1βˆ’πœ–2ξ€Έξ€·(1βˆ’πœ‡)1βˆ’π‘„πœ‡πœ–2ξ€Έξ€·1βˆ’πœ‡πœ–2ξ€Έ2,(5.19) where we introduced πœ–=π‘’βˆ’π‘‡/21,πœ‡=π‘„π‘’πœ.(5.20) Substituting these end-points back to the formula for the resolvent, we find πœ”Β±1(𝑝)=𝑝𝑇log1+πœ‡πœ–2𝑝+1+π‘„πœ–2ξ€Έβˆ“ξ€·1βˆ’πœ‡πœ–2ξ€Έβˆš(π‘βˆ’π‘Ž)(π‘βˆ’π‘)2π‘’βˆ’π‘‡(ξƒͺ𝑝+𝑄).(5.21) As a check, this result indeed satisfies the consistency condition (3.20) πœ”+(𝑝)+πœ”βˆ’πœ•(𝑝)=π‘π‘‰πœ(𝑝)𝑇,(5.22) with π‘‰πœ given in (5.16). From the knowledge of the resolvent, we can also determine eigenvalue density along the cut (3.21) 1𝜌(𝑝)=𝑝𝑇log1+πœ‡πœ–2𝑝+1+π‘„πœ–2βˆ’ξ€·1βˆ’πœ‡πœ–2ξ€Έβˆš(π‘βˆ’π‘Ž)(π‘βˆ’π‘)ξ€·1+πœ‡πœ–2𝑝+1+π‘„πœ–2+ξ€·1βˆ’πœ‡πœ–2ξ€Έβˆšξƒͺ(π‘βˆ’π‘Ž)(π‘βˆ’π‘),(5.23) as well as the spectral curve. Writing π‘₯=π‘π‘‡πœ”(𝑝) and 𝑝=𝑒𝑦, and after a few simple rescalings, we find that the spectral curve takes form 𝑒π‘₯+𝑦+𝑒π‘₯+𝑒𝑦+𝑄1𝑒2π‘₯+𝑄2𝑒2𝑦+𝑄3=0,(5.24) where 𝑄1=πœ–21+πœ‡π‘„ξ€·1+πœ‡πœ–2ξ€Έξ€·1+π‘„πœ–2ξ€Έ,𝑄2=πœ‡1+π‘„πœ–2ξ€·(1+πœ‡π‘„)1+πœ‡πœ–2ξ€Έ,𝑄3=𝑄1+πœ‡πœ–2ξ€·1+πœ–2𝑄.(1+πœ‡π‘„)(5.25)

The above curve is given by a symmetric function of 𝑄, πœ‡=π‘„βˆ’1π‘žπ‘› and πœ–2=π‘’βˆ’π‘‡. Apparently this is a mirror curve of the so-called closed topological vertex geometry, which is a Calabi-Yau manifold arising from a symmetric resolution of β„‚3/β„€2Γ—β„€2 orbifold, see Figure 15. This geometry has three moduli, that is, the original KΓ€hler moduli 𝑄 of the resolved conifold, the chamber parameter 𝑛 (encoded in πœ‡), and the 't Hooft parameter 𝑇, which are all unified in a geometric way in our matrix model. Moreover, the fractional coefficients in the curve equation (5.24) encode the correct mirror map for the closed topological vertex geometry (and to the linear order, these coefficients are just 𝑄, πœ‡, and π‘’βˆ’π‘‡).

In the BPS counting problem we are interested in, as follows from the form of the identity operator (5.1), ultimately we need to consider matrices of infinite size. We also need to keep fixed 𝑔𝑠, so we should consider the limit of π‘‡β†’βˆž, or equivalently πœ–β†’0. Up to a linear shift of π‘₯ and 𝑦, the (5.24) in this limit becomes πœ‡π‘’2𝑦+𝑒π‘₯+𝑦+𝑒π‘₯+(1+π‘„πœ‡)𝑒𝑦+𝑄=0.(5.26) The manifold corresponding to this curve is the suspended pinch point (SPP) geometry, with 𝑄 and πœ‡ encoding flat coordinates representing sizes of its two β„™1's. Having found the SPP mirror curve, let us also make the following remarks.

Firstly, we see that not only the spectral and mirror curves agree, but moreover the matrix integral (5.9) reproduces (after the identification