Vortex core order and field-driven phase coexistence in the attractive Hubbard model

Vortex core order and field-driven phase coexistence in the attractive Hubbard model

Madhuparna Karmakar madhuparna.k@gmail.com The Institute of Mathematical Sciences, HBNI, C I T Campus, Chennai 600 113, India    Gautam I. Menon menon@imsc.res.in The Institute of Mathematical Sciences, HBNI, C I T Campus, Chennai 600 113, India    R. Ganesh ganesh@imsc.res.in The Institute of Mathematical Sciences, HBNI, C I T Campus, Chennai 600 113, India
July 18, 2019

Superconductivity occurs in the proximity of other competing orders in a wide variety of materials. Such competing phases may reveal themselves when superconductivity is locally suppressed by a magnetic field in the core of a vortex. We explore the competition between superconductivity and charge density wave order in the attractive Hubbard model on a square lattice. Using Bogoliubov-deGennes mean field theory, we study how vortex structures form and evolve as the magnetic flux is tuned. Each vortex seeds a CDW region whose size is determined by the energy cost of the competing phase. The vortices form a lattice whose lattice parameter shrinks with increasing flux. Eventually, their charge-ordered vortex cores overlap, leading to a field-driven coexistence phase exhibiting both macroscopic charge order and superconductivity – a ‘supersolid’. Ultimately, superconductivity disappears via a first-order phase transition into a purely charge ordered state. We construct a phase diagram containing these multiple ordered states, using , the next-nearest neighbour hopping, to tune the competition between phases.



Superconductivity is often obtained in proximity to other ordered ground states. The most prominent example being the high T cuprates, where superconductivity competes with antiferromagnetism and with charge orderLake et al. (2001); Chang et al. (2012a). A particularly interesting way to stabilize underlying competing phases is to apply a magnetic field, locally suppressing superconductivity to create vortices. The core region of the vortex can then host competing correlationsHoffman et al. (2002); Curran et al. (2011); Machida et al. (2016). Indeed, experiments with scanning tunnelling microscopy have revealed charge-orderedHoffman et al. (2002); Machida et al. (2016) vortex cores in the cuprates. NMR studies of YBaCuO indicate that as the magnetic field increases, the inter-vortex distance decreases; at a critical field strength, vortex cores overlap leading to charge order throughout the systemWu et al. (2013). These and related experiments motivate the study of vortex core order and field-driven coexistence in the attractive Hubbard model, the simplest model to show competition between superconductivity (SC) and charge density wave (CDW) order.

Hubbard model and symmetry:

We consider fermions on a square lattice, described by


where is the chemical potential and is the strength of the on-site attractive interaction (). The hopping parameter takes the value (henceforth set to unity) for nearest neighbours and is zero otherwise. When is tuned to half-filling, this model possesses a remarkable symmetry with SC and CDW order becoming degenerate; the order parameters form an enlarged space having symmetry as shown in Fig. 1(left)Yang (1989); Zhang (1990); Yang and Zhang (1990); Burkov and Paramekanti (2008); Ganesh et al. (2009). This is a delicate symmetry arising from the bipartite nature of the square lattice with hoppings connecting sites of different sublattices. We can tune away from this degenerate point by introducing a next-nearest neighbour hopping, . The term lowers the energy of the SC phase relative to CDW phase.

Figure 1: Left: The space of order parameters forming an sphere; the equator corresponds to the phase of the SC order parameter while the poles correspond to two possible checkerboard CDW orders. A generic point on the sphere represents coexisting SC and CDW orders. Right: The order parameters forming a ‘meron’ in the vicinity of a vortex, with . Far from the core, the pseudospins lie in the plane and wind by as we move around the vortex. Within the core, they cant out of the plane to generate CDW order.

The degeneracy leads to a local pseudospin order parameter whose components are as shown in Fig. 1. Here, and (defined below) are the local superconducting and CDW order parameters. This symmetry is directly analogous to the hypothesized symmetryDemler et al. (2004) in the cuprates which groups SC and antiferromagnetism into an enlarged order parameter space. As a testable consequence of theory, it was proposed that vortex cores would have antiferromagnetic orderArovas et al. (1997); Hu and Zhang (2002). Analogously, the Hubbard model in Eq. 1 will possess CDW order in the vortex core. In the language of pseudospins, a vortex corresponds to a ‘meron’, as shown in Fig. 1(right) – in the core region, the moments cant out of the plane to locally give rise to CDW order. Here, unlike in the cuprates, we have direct control over the symmetry breaking in the form of the hopping. We study the Hubbard model in an applied field demonstrating CDW order in the vortex core. Our key result is a field-driven SC-CDW coexistence regime which arises from the overlap of vortex cores. This state simultaneously breaks translational symmetry and -gauge symmetry, demonstrating a new route to ‘supersolidity’Boninsegni and Prokof’ev (2012).

Figure 2: (a) Superconducting and (b) CDW order profiles at different . The inset to panel (a) shows the underlying length scale vs. . The inset to panel (b) shows the FWHM widths, and , vs. . The lower panels show the spatial maps of SC() (left) and density()(right) order parameters around a vortex core for . The interaction strength is fixed at .

Bogoliubov deGennes mean-field theory:

We perform simulations on an lattice with periodic boundary conditions, with up to . To introduce an orbital magnetic field, we add a complex phase to the hopping amplitudes given by , where is the vector potential, see Supplementary Materials for details. The net magnetic flux through a closed surface must be quantized in units of Dirac (1931). We take the net flux through our system to be where is an integer. As each vortex carries a flux , we will always have an even number of vortices in the system. In particular, the lowest magnetic flux we can have is , corresponding to two vortices.

We decompose the on-site interaction term in pairing and density channels. The SC order parameter is complex-valued, defined as . The density order parameter is defined as . The local CDW order parameter can be defined as , which measures the local deviation from half-filling. With these mean-field parameters, the Hamiltonian takes the form of a matrix, which can be diagonalized using the Bogoliubov-Valatin transformation Chen and Ting (2002, 2005); Han (2010). We obtain self-consistent values of and on every site.

We find several self-consistent mean-field configurations, of which the one with lowest energy is to be chosen. One solution is a pure CDW state in which for all and , corresponding to uniform CDW order. In the absence of a magnetic field and in the presence of a non-zero , this state has higher energy than the uniform SC phase. When a field is imposed, this state is not affected as it is insulating – its energy remains constant, independent of the flux (see Supplementary Materials). In contrast, the SC phase necessarily develops vortices when a field is imposed. As the number of vortices increases with flux, so does the energy of the SC. As seen from these energetic arguments, an applied magnetic field induces competition between SC and CDW orders.

The symmetry of the attractive Hubbard model only exists precisely at half-filling. As we are interested in phase competition, all results presented here are at half-filling. We present results for for the following reason. At large , the Hubbard model can be mapped to a spin problem with antiferromagnetic superexchange interactionsBurkov and Paramekanti (2008); Ganesh et al. (2009). The local order parameter is, in fact, the spin whose components are as shown in Fig. 1. At low temperatures, we expect the system to have uniform spin length, ie., , a constant independent of position. SC and CDW order parameters are not independent, having to satisfy this uniform-length constraintArovas et al. (1997). With these considerations, the appropriate Landau Ginzburg free energy density is given byArovas et al. (1997)


The order parameters are coupled by the uniform length constraint: . SC and CDW orders become degenerate when and the magnetic field is turned off, revealing the underlying symmetry. At , we find that the mean-field results always satisfy the uniform spin length constraint and the above Landau Ginzburg theory applies. We find the same qualitative results extending to small values as well.

Motivated by recent experiments revealing charge order in the cuprates, several authors have studied field theories similar to Eq. 2Efetov et al. (2013); Hayward et al. (2014); Wachtel and Orgad (2014) with Ref. Meier et al., 2013 also incorporating an orbital magnetic field. Our study of the attractive Hubbard model at strong coupling can be viewed as an ultraviolet regularization of such a field theory.

Vortex profile:

Setting , we obtain the lowest flux configuration with two well-separated vortices. As is increased, we find CDW order in the vortex core until . For larger values, we find a normal core with no CDW correlations. Figs. 2(a) and 2(b) show the profiles of superconducting and CDW order at selected values of . The lower panels of Fig. 2 show the spatial maps of SC and CDW order parameters around a single vortex for : CDW correlations can be clearly seen in the vortex core region. The same information is presented in spin language in Fig. 1(right).

The SC and CDW profiles are, in fact, set by the same length scale, , as the order parameters satisfy the uniform spin-length constraint. We obtain by fitting the SC profile to ; this functional form is consistent with the free energy in Eq. 2Arovas et al. (1997). The resulting is plotted in the inset to Fig. 2(a). Separately, we define two length scales, and , as the full-widths at half-maximum of SC and CDW profiles respectively. We find that is always larger than as shown in the inset to Fig. 2(b), although there is only one underlying length scale, . Superconductivity is suppressed in the vortex core out to a radius set by . However, each vortex also hosts CDW correlations which extend to a larger distance, . This sets the stage for a coexistence phase, which can be understood as follows. As the field is increased, the vortices are packed more and more tightly. Naïvely, superconductivity persists until the inter-vortex distance approaches . Much before this, when the inter-vortex distance reaches , the CDW regions around each vortex overlap. CDW order percolates throughout the system on top of the SC background, leading to a field-driven coexistence phase as we demonstrate below. Similar arguments have recently been put forward in the cuprates based on NMR resultsWu et al. (2013).

Figure 3: Spatial maps of the SC order parameter amplitude, Fourier transform of SC amplitude and the density order parameter. The results are for , and five different magnetic field values, parametrized by . We also show the density of states of fermionic excitations, , as a function of field.

Vortex lattice evolution:

Fig. 3 shows our results for with varying (the total magnetic flux being ) on a lattice. The panels show real space maps of and , showing the evolution of a vortex lattice with increasing flux. In addition, we plot the Fourier transform of the SC order parameter, defined as . The distribution of peaks in reveals the geometry of the vortex lattice. The figure also plots the electronic density of states in the mean field ground state.

For small fields, with , we find well separated vortices forming an anisotropic triangular lattice. At , the lattice becomes near-isotropic. Upon increasing the field to , the vortices form a square lattice. This suggests an underlying phase transition driven by tuning vortex density. A similar transition into a second vortex lattice phase has been suggested in YBCO from torque magnetometry resultsYu et al. (2016). For , we find phase separation into square and triangular vortex lattices. Finally, at , we find a first order phase transition to a pure CDW phase, which has lower energy than solutions with SC order. Thus, at mean-field level, is set by the competing CDW phase, unlike in conventional superconductors.

Figure 4: CDW parameter (red diamonds) and the quasiparticle gap (black circles) as a function of flux at as obtained from mean-field theory. The CDW regions around vortex cores begin to overlap at , heralding coexistence. The gap starts to increase beyond this threshold. In the high field ‘CDW’ region, a pure CDW mean field state is favoured over a SC phase.

Phase coexistence:

As seen in Fig. 3, every vortex core nucleates CDW correlations which begin to overlap when . The CDW order becomes progressively stronger with increasing field as vortex cores overlap more and more. For , we have near-uniform CDW order. With increasing field, the SC order weakens while the CDW order parameter grows. As a result, the electronic gap never closes, as shown in Fig. 3.

Fig. 4 shows the in-field phase diagram for . As we force a magnetic flux through Peierl’s substitution, there is no in our simulations. To quantify the strength of CDW order, we define the parameter where is the Fourier component of the density order parameter: . In a pure CDW state with site occupation oscillating between 0 and 2, takes the value unity. For small values, we find to be small, indicating weak CDW order arising from well separated vortex cores. With increasing , increases monotonically. Within our mean-field theory, we find that CDW correlations begin to span the system when . Based on this observation, we use as a heuristic criterion to signal macroscopic CDW order.

While CDW and SC compete spatially, they both serve to open an electronic gap. At zero magnetic field, we have a uniform SC state with a large gap of the order of . With an applied field, we induce vortices with CDW correlations, leading to a spatially textured order parameter field. Initially, for small magnetic fields, the order parameter gradients reduce the electronic gap. Once CDW order percolates throughout the system, the CDW order parameter no longer suffers sharp gradients and strengthens the gap once again.


Working in the context of the attractive Hubbard model, we have shown that competing CDW order emerges in vortex cores. At large fields, vortex cores overlap leading to a coexistence phase with SC and CDW correlations spanning the system – a lattice version of a supersolid. The existence of a supersolid has been heavily debated in the context of liquid HeBoninsegni and Prokof’ev (2012). Our study suggests that superconductors with competing phases are strong candidates for supersolidity.

Our mean field theory in the strong coupling limit can also be seen as a spin problem with chiral interactions introduced by the orbital magnetic field. The vortex lattice ground state can be interpreted as a ‘meron crystal’ – a pseudospin state with a chiral textureVolovik and Krusius (2012). When CDW order spans the system, this state spontaneously breaks a symmetry corresponding to two possible checkerboard density patterns, or equivalently to the choice between and ordering of the spins. Our estimate for the superfluid stiffness (see Supplementary Material) indicates that this state is stable to fluctuations at intermediate fields. At large magnetic fields, the CDW order becomes much stronger than SC. In this regime, fluctuations may destabilize SC while leaving the CDW order intact. This suggests a ‘vortex liquid’ phase with remnant CDW order and vanishing superfluid stiffness. This is an exciting direction for future study. In fact, field-induced coexistence and a charge-ordered vortex-liquid state have been reported in YBCOLeBoeuf et al. (2013); Gerber et al. (2015); Yu et al. (2016).

We present a schematic phase diagram for the Hubbard model in the - plane in Fig. 5. The plot shows the region, where is the critical value beyond which CDW order is energetically unfavourable at all magnetic fields. This phase diagram provides a reference point for understanding phase competition in the cuprates. For example, we find that SC is lost at via a first order transition to the competing CDW phase; this may have been seen in YBCO using thermal conductivity measurementsGrissonnanche et al. (2014). We see that vortices themselves show interesting ordering phenomena with the vortex lattice changing from triangular to square configuration. This resembles the suggestion of two vortex solid states in YBCO from torque magnetometry measurementsYu et al. (2016).

Figure 5: Schematic phase diagram at in the plane. The lower phase boundary represents a crossover field at which CDW correlations begin to span the system. The upper boundary represents , a first order phase transition into a pure CDW state.

While our study focusses on a simple model Hamiltonian, our results broadly apply to several material families which host competing orders. We have used as a convenient handle to tune phase competition, this role could be played by experimentally tunable parameters such as doping in the cuprates Chang et al. (2012b), pressure in TiSeKusmartseva et al. (2009), etc. Our results provide a theoretical paradigm to understand phase competition in these systems.
Acknowledgments: We thank Arun Paramekanti, David Hawthorn, David Broun and Anton Burkov for discussions. MK acknowledges the use of High performance computing facility at Harish Chandra Research Institute, Allahabad, India.


Appendix A Supplementary Material

a.1 Peierl’s substitution

We have two types of hopping terms in our Hamiltonian: nearest-neighour hopping with amplitude and next-nearest hopping with amplitude . The phase of each hopping element represents a line integral of the vector potential, in accordance with the principle of Peierl’s substitution. As the vector potential is not uniquely defined, there are several possible ways to assign the complex phases. The gauge invariant quantity is the magnetic flux: the sum of the hopping-phases along closed loops on the lattice.

In our periodic lattice, we assign hopping-phases so as to obtain a uniform magnetic flux. We use the scheme shown in Fig. S1. We have introduced a parameter which encodes the phase picked up by an electron when hopping around any square plaquette, i.e., the magnetic flux through each square plaquette is . Our square lattice system with periodic boundary conditions is equivalent to a torus, a closed surface. As shown by Dirac, the magnetic flux through any closed surface must be quantized in units of so that , with being an integer. The parameter determines the total flux through the system, given by .

Dirac further argued that in the presence of a non-zero flux, we cannot define electronic wavefunctions on the surface smoothly. The phase of the wavefunctions must wind around a singularity, which is called the ‘Dirac string’. In our scheme, the Dirac string passes through the red region within the top-right square plaquette in Fig. S1. The sum of the hopping-phases around a contour which encloses this region has an additional contribution of . This does not indicate an increased magnetic flux through this region; rather, it reflects a singularity in the definition of the electronic wavefunctions. Indeed, as this phase is a multiple of , it does not lead to any observable consequences.

Figure S1: Peierls substitution scheme on the lattice. Left: The periodic cluster for for illustration purposes. Sites are labelled as – not all site labels are shown. The bonds present in the cluster are depicted using dark solid lines. Bonds repeated due to periodic boundary conditions are shown in dotted lines. The square plaquettes are divided into two types as indicated by the colour. Right: Representative plaquettes of each type are shown, with the lower-left site labelled. The hopping-phases for the bonds on these plaquettes are shown. The parameter determines the flux through each square plaquette, given by . The Dirac string passes through the red triangular area. Due to the Dirac string, the flux is constrained to satisfy , where is an integer.

a.2 Bogoliubov-de-Gennes formalism

Figure S2: Phase competition between SC and CDW orders. (Left) Comparison of ground state energies as a function of , with the magnetic field turned off (). At , the two orders are degenerate. A non-zero value of lowers the energy of the SC state. (Right) Ground state energies as a function of magnetic flux, , with fixed at . At , the CDW state wins over the SC state indicating a first order phase transition. All energies are calculated at half-filling.

The Hamiltonian for the Hubbard model is given in the main text. We decompose the on-site interaction term in pairing and density channels via a mean field decomposition. The complex superconducting order parameter is defined as, , while the charge order parameter is defined as, . The resulting effective Hamiltonian is given by


The hopping phases and are assigned according to the Peierl’s substitution scheme described above. We diagonalize this Hamiltonian using a Bogoliubov-Valatin transformation given by , where () creates (annihilates) a quasiparticle with spin with energy and wavefunctions and . We have introduced a spin index and . The resulting gap and number equations are

where if is the Fermi function at zero temperature. Starting from initial guess values, we iterate these equations to obtain self-consistent values of and . The chemical potential is tuned to fix the density at half-filling.

Figure S3: Spatial maps of the pairing amplitude () and the CDW () order with changing magnetic field at on a 2424 lattice.

a.3 Phase competition

For a given set of parameters , and , we find several self-consistent solutions. In particular, we find a pure CDW state with . To illustrate the phase competition in the Hubbard model, we compare the energy of this CDW state with that of the SC state in Fig. S2. In the absence of a magnetic field (), the two states are degenerate when , while a non-zero lowers the energy of the SC state.

Competition between orders can also be tuned by increasing the magnetic field. This is depicted in Fig. S2(right). When is increased at fixed , the energy of the CDW state does not change as the CDW state is an insulator. However, in the SC phase, increasing introduces more vortices and increases the energy. The energy of the SC state steadily rises and eventually crosses the CDW energy. This signals at the mean-field level, with CDW order becoming energetically favourable over a superconducting vortex state.

a.4 Higher regime

With increasing , the radius of the CDW region in each vortex core shrinks. Consequently, the threshold field for coexistence increases. Fig. S3 shows spatial maps of the SC and CDW order parameters for which can be compared with the data in Fig. 3 of the main text. The percolation of CDW correlations is slow to occur with a coexistence state only setting at . Finally, is encountered at , when the CDW state becomes energetically favourable. Unlike the case of discussed in the main text, it is difficult to discern changes in the geometry of the vortex lattice here due to the high density of vortices.

a.5 Superfluid stiffness

Figure S4: Superfluid stiffness. We plot vs. for three sets of values: (a) (blue downward triangles), (b) (green upward triangles), and (c) (red squares). The lines are fits to the form . All three best-fit lines have positive slopes indicating a positive superfluid stiffness.

Our mean-field results indicate coexistence of SC and CDW orders, forming a supersolid state. To check if this phase is stable to fluctuations, the standard diagnostic is superfluid stiffness – which measures the energy cost of imposing a smooth gradient in the SC order parameter. To estimate the stiffness, we take the following route. We introduce an additional component of the vector potential . If we were to turn off the orbital magnetic field, this vector potential leads to a ‘flowing’ superfluid solution with . The resulting energy cost is a measure of superfluid stiffness. For a simple superfluid with no competing order, this energy cost (the increase in energy per site) is proportional to , where is the superfluid stiffness and is the system size. We define , where is the energy (per site) of the flowing state. This is calculated by using the values obtained after inclusion of the tangential vector potential to evaluate the expectation value of the mean-field Hamiltonian.

The obtained values are plotted as a function of in Fig. S4. In the case (no orbital magnetic field), we see that indeed scales as , with a positive slope. This slope is proportional to the superfluid stiffness.

In the presence of the orbital field, we seek to plot for configurations with the same magnetic flux density across different system sizes. In our calculations, the magnetic flux density is , with being an integer and ranging from (for smaller sizes, we see strong finite size effects). Generically, it is not possible to find multiple values for which is a constant. In Fig. S4, we plot for and which correspond to approximately constant values. The resulting values also scale linearly with with a positive slope. We conclude that these superfluid stiffness is positive for these flux densities. We also note that the stiffness decreases with increasing flux density. In particular, we note that the stiffness is positive for the flux densities corresponding to and on a 2424 lattice. As discussed in the main text, these parameters have macroscopic CDW order with CDW correlations spanning the entire system. This suggests that the coexistence phase is stable to fluctuations.

It is possible that that stiffness may vanish at higher flux densities, perhaps close to . In this regime, the CDW order parameter becomes approximately constant while the SC order parameter suffers large gradients due to the presence of vortices. It is then conceivable that fluctuations can wash out the in-plane order while preserving order in the -direction. This will lead to a ‘pairing liquid’ state (analogous to a spin liquid) with remnant CDW ordering. Equivalently, this can be understood as melting of the vortex lattice. Below this threshold, the inter-vortex distances are fixed by strong interactions between vortices. A small amount of disorder will then suffice to pin the entire vortex lattice so as to generate a robust SC state. However, when fluctuations wash out the coherence in the SC, the vortices become mobile giving rise to a ‘vortex liquid’. This is an interesting direction for future study.

Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description