Disordered ultracold atomic gases in optical lattices: A case study of Fermi-Bose mixtures

V. Ahufinger, L. Sanchez-Palencia, A. Kantian, A. Sanpera,* and M. Lewenstein* Grup d’Òptica, Departament de Física, Universitat Autònoma de Barcelona, E-08193 Belaterra, Barcelona, Spain Institut für Theoretische Physik, Universität Hannover, D-30167 Hannover, Germany Laboratoire Charles Fabry, Institut d’Optique Théorique et Appliquée, Université Paris-Sud XI, F-91403 Orsay Cedex, France Institut für Quantenoptik und Quanteninformation der Österreichischen, Akademie der Wissenschaften, A-6020 Innsbruck, Austria Institut für Theoretische Physik, Universität Innsbruck, A-6020 Innsbruck, Austria Grup de Física Teòrica, Departament de Física, Universitat Autònoma de Barcelona, E-08193 Belaterra, Barcelona, Spain Institut de Ciències Fotòniques, E-08034 Barcelona, Spain Received 7 September 2005; published 21 December 2005

Since the discovery of the quantum localization phenomenon by Anderson in 1958 ͓1͔, disordered and frustrated systems have played a central role in condensed matter physics. They have been involved in some of the most challenging open questions concerning many body systems ͑cf. ͓2-4͔͒. Quenched ͑i.e., frozen on the typical time scale of the considered systems͒ disorder determines the physics of various phenomena, from transport and conductivity, through localization effects and metal-insulator transition ͑cf. ͓5͔͒, to spin glasses ͑cf. ͓6-8͔͒, neural networks ͑cf. ͓9͔͒, percolation ͓10,11͔, high T c superconductivity ͑cf. ͓12͔͒, or quantum chaos ͓13͔. One of the reasons why disordered systems are very hard to describe and simulate is related to the fact that, usually, in order to characterize the system, one should calculate the relevant physical quantities averaged over a particular realization of the disorder. Analytical approaches require the averaging of, for instance, the free energy, which ͑being proportional to the logarithm of the partition function in the canonical ensemble͒ is a very highly nonlinear function of the disorder. Averaging requires then the use of special methods, such as the replica trick ͑cf. ͓6͔͒, or supersymmetry method ͓14͔. In numerical approaches this demands either simulating very large samples to achieve "selfaveraging," or numerous repetitions of simulations of small samples. Obviously, this difficulty is particularly important for quantum disordered systems. Systems which are not disordered but frustrated ͑i.e., unable to fulfill simultaneously all the constrains imposed by the Hamiltonian͒, lead very often to similar difficulties, because quite often they are char-acterized at low temperature by an enormously large number of low energy excitations ͑cf. ͓15͔͒. It is thus desirable to ask whether atomic, molecular physics, and quantum optics may help to understand such systems. In fact, very recently, it has been proposed how to overcome the difficulty of quenched averaging by encoding quantum mechanically in a superposition state of an auxiliary system, all possible realizations of the set of random parameters ͓16͔. In this paper we propose a more direct approach to the study of disorder: direct realization of various disordered models using cold atoms in optical lattices.

B. Disordered ultracold atomic gases
In recent years there has been enormous progress in the studies of ultracold weakly interacting ͓17͔, as well as strongly correlated, atomic gases. In fact, present experimental techniques allow one to design, realize, and control in the laboratory various types of ultracold interacting Bose or Fermi gases, as well as their mixtures. Such ultracold gases can be transferred to optical lattices and form a, practically perfect, realization of various Hubbard models ͓18͔. This observation, suggested in the seminal theory paper by Jaksch et al. ͓19͔ in 1998, and confirmed then by the seminal experiments of Greiner et al. ͓20͔, has triggered enormous interest in the studies of strongly correlated quantum systems in the context of atomic and molecular lattice gases.
It soon became clear that one can introduce local disorder and/or frustration to such systems in a controlled way using various experimentally feasible methods. Local quasidisorder potentials may be created by superimposing superlattices incommensurable with the main one to the system. Although strictly speaking such a superlattice is not disordered, its effects are very similar to those induced by the genuine random potentials ͓21-23͔. Controlled local truly random poten-*Also at Institució Catalana de Recerca i Estudis Avançats. tials can be created by placing a speckle pattern on the main lattice ͓24,25͔. As shown in Ref. ͓21͔, for a system of strongly correlated bosons located in such a disordered lattice, both methods should permit one to achieve an Anderson-Bose glass ͓26͔, provided that the correlation length of the disorder L dis is much smaller than the size of the system. Unfortunately, it is difficult to have L dis smaller than a few microns using speckles. Thus, L dis is typically larger than the condensate healing length l heal =1/ ͱ 8na, where n is the condensate density, and a the atomic scattering length. Due to this fact, i.e., due to the effects resulting from the nonlinear interactions, it is difficult to achieve the Anderson localization regime with weakly interacting Bose-Einstein condensates ͑BECs͒ ͓27,28͔. We have shown, however, that quantum localization should be experimentally feasible using the quasidisorder created by several lasers with incommensurable frequencies ͓29͔. Random local disorder appears also, naturally, in magnetic microtraps and atom chips as a result of roughness of the underlying surface ͓͑30͔, for theory see ͓31,32͔͒.
One can also create disorder using a second atomic species, by rapidly quenching it from the superfluid to the localized Mott phase. After such a process, different lattice sites are populated by a random number of atoms of the second species, which act effectively as random scatterers for the atoms of the first species ͓33͔. Last, but not least it is possible to use Feschbach resonances in fluctuating or inhomogeneous magnetic fields in order to induce a type of disorder that corresponds to random, or at least inhomogeneous, nonlinear interaction couplings ͓34͔ ͑for theory in onedimensional ͑1D͒ systems see ͓35,36͔͒. It has been also been proposed ͓37͔ that tunneling induced interactions in systems with local disorder results in controllable disorder on the level of next-neighbor interactions. That opens a possible path for the realization of quantum spin glasses ͓37͔. As we have already mentioned, several experimental groups have already achieved ͓27-29͔, or are soon going to realize ͓33,34͔ disordered potentials using these methods. It is worth mentioning here a very recent attempt to create controlled disorder using optical tweezers methods ͓38͔. There are also several ways to realize nondisordered but frustrated systems with atomic lattice gases. One is to create such gases in "exotic" lattices, such as the Kagomé lattice ͓39͔, another is to induce and control the nature and range of interactions by adjusting the external optical potentials, such as, for instance, proposed in Ref. ͓40͔. Another example of such a situation is provided for instance by atomic gases in a two-dimensional lattice interacting via dipole-dipole interactions with dipole moment polarized parallel to the lattice ͓41͔.
Finally, there are also several ideas concerning the possibility of realizing various types of complex systems using atomic lattice gases or trapped ion chains. Particularly interesting here are the possibilities of producing long range interactions ͑falling off as inverse of the square, or cube of the distance͒ ͓42͔, systems with several metastable energy minima, and last, but not least systems in designed external magnetic ͓43͔, or even non-Abelian gauge fields ͓44͔.

C. Quantum information with disordered systems
One important theoretical aspect that should be considered in this context deals with the role of entanglement in quantum statistical physics in general ͑where it concerns quantum phase transitions, entanglement correlation length, and scaling ͓45͔͒, and characterization of various types of distributed entanglement. This aspect is particularly interesting in theoretical and experimental studies of disordered systems. The question which one is tempted to ask is, to what extent can one realize quantum information processing in ͑i͒ disordered systems, ͑ii͒ nondisordered systems with long range interactions, and ͑iii͒ nondisordered frustrated systems.
At the first sight, the answer to this question is negative. Disordered quantum information processing sounds like contradictio in adjecto. But, one should not neglect the possible advantages offered by the systems under investigation. First, such systems have typically a significant number of ͑local͒ energy minima, as, for instance, happens in spin glasses. Such metastable states might be employed to store information distributed over the whole system, as in neural network models. The distributed storage implies redundancy similar to the one used in error correcting protocols ͓46͔. Second, in the systems with long range interactions the stored information is usually robust: metastable states have large basins of attraction thermodynamically, and destruction of a part of the system does not destroy the metastable states ͑for the preliminary studies see Refs. ͓47,48͔͒. Third, and perhaps the most important aspect for the present paper, is that atomic ultracold gases offer a unique opportunity to realize special purpose quantum computers ͑quantum simulators͒ to simulate quantum disordered systems. The importance of the experimental realizations of such quantum simulators will without doubts forward our understanding of quantum disordered systems enormously. In particular, we can think about large scale quantum simulations of the Hubbard model for spin 1 / 2 fermions with disorder, which lies at the heart of the present-day-understanding of high T c superconductivity. The impact of this possibility for physics and technology is hard to overestimate. Fourth, very recently, several authors have used the ideas of quantum information theory to develop novel algorithms for efficient simulations of quantum systems on classical computers ͓49͔. Applications of these novel algorithms to disordered systems are highly desired.

D. Fermi-Bose mixtures
The present paper deals with the above formulated questions, which lie at the frontiers of modern theoretical physics, and concern not only atomic, molecular, and optical ͑AMO͒ physics and quantum optics, but also condensed matter physics, quantum field theory, quantum statistical physics, and quantum information. This interdisciplinary theme is one of the most hot current subjects of the physics of ultracold gases. In particular, we present here the study of a specific example of disordered ultracold atomic gases: Fermi-Bose ͑FB͒ mixtures in optical lattices in the presence of additional inhomogeneous and random potentials.
In the absence of disorder and in the limit of strong atomatom interactions such systems can be described in terms of composite fermions consisting of a bare fermion, or a fermion paired with one boson ͑bosonic hole͒, or two bosons ͑bosonic holes͒, etc. ͓50͔. The physics of Fermi-Bose mixtures in this regime has been studied by us recently in a series of papers ͓52-54͔; for contributions of other groups to the studies of FB mixtures in traps and in optical lattices see Ref.
͓55͔ and for the studies of strongly correlated FB mixtures in lattices see Ref. ͓63͔, respectively. In particular, the validity of the effective Hamiltonian for fermionic composites in 1D was studied using exact diagonalization and the density matrix renormalization group method in Ref. ͓64͔. The effects of inhomogeneous trapping potential on FB lattice mixtures has been for the first time discussed by Cramer et al. ͓65͔. The physics of disordered FB lattice mixtures was studied by us in Ref. ͓37͔, which has essentially demonstrated that this system may serve as a paradigm fermionic system to study a variety of disordered phases and phenomena: from Fermi glass to quantum spin glass and quantum percolation.

E. Plan of the paper
The main goal of the present paper is to present the physics of the disordered FB lattice gas in more detail, and in particular to investigate conditions for obtaining various quantum phases and quantum states of interest.
The paper is organized as follows. Section II describes the "zoology" of disordered systems and disordered phases known from condensed matter physics. We pay particular attention to the systems realizable with cold atoms on one side, and particularly interesting from the other. This latter phrase means that we consider here the systems that concern important open questions and challenges of the physics of disordered systems. In this sense this section is thought as a list of such challenging open questions that can be perhaps addressed by the cold gases community. This section is thus directed to the cold gases experts, and is supposed to motivate and stimulate their interest in the physics of ultracold disordered systems.
In Sec. III, we introduce the composite fermions formalism, first discussing it for the case of homogeneous lattices, and then for disordered ones. We derive here the explicit formulas for the effective Hamiltonian, and for various types of disorder. One of the results of this section concerns the generalizations of the results of Ref. ͓37͔, that implies that local disorder on the level of the Fermi-Bose Hubbard model leads to randomness of the nearest neighbor tunneling and coupling coefficients for the composite fermions. Obviously, these tunneling and coupling coefficients arise in effect of tunneling mediated interactions between the composites.
In Sec. IV, we present our numerical results in the weak disorder limit, based on the time dependent Gutzwiller ansatz. These results concern mainly the physics of composites, the realization of Fermi glass, and the transition from Fermi liquid to Fermi glass.
The results for the case of strong disorder, spin glasses, are discussed in Sec. V. The problem of detection of the phenomena predicted in this paper is addressed in Sec. VI, whereas we conclude in Sec. VII.

II. DISORDERED SYSTEMS: THE OLD AND NEW CHALLENGES FOR AMO PHYSICS
In this section we present a list of problems and challenges of the physics of disordered systems that may, in our opinion, be realized and addressed in the context of physics of ultracold atomic or molecular gases. We concentrate here mainly on fermionic systems. This section is written on an elementary level and addressed to nonexperts in the physics of disordered systems.

A. Anderson localization
One of the most spectacular effects of disorder concerns single particles. The spectrum of a Hamiltonian of a free particle in free space or in a periodic lattice is continuous and the corresponding eigenfunctions are extended ͑plane waves or Bloch functions͒. Introduction of disorder may drastically change this situation. The basic knowledge about these phenomena comes from the famous scaling theory formulated by the "gang of four" ͓͑66,5͔͒.
The scaling theory predicts that in 1D infinitesimally small disorder leads to exponential localization of all eigenfunctions. The localization length ͑defined as the width of the exponentially localized states͒ is a function of the ratio between the potential and the kinetic ͑tunneling͒ energies of the eigenstate and the disorder strength. For the case of discrete systems with constant tunneling rates and local disorder distributed according to a Lorentzian distribution ͑Lloyd's model, cf. ͓13͔͒ the exact expression for the localization length is known. In general, an exact relation between the density of states and the range of localization in 1D has been provided by Thouless ͓67͔. Hard core bosons with strong repulsion in 1D chains, described by the XY model in a random transverse field, can be mapped using the Wigner-Jordan transformation to 1D noninteracting fermions in a random local potential, which in turn maps the bosonic problem onto the problem of Anderson localization ͓68͔.
In 2D, following the scaling theory, it is believed that localization occurs also for arbitrarily small disorder, but its character interpolates smoothly between algebraic for weak disorder, and exponential for strong disorder. There are, however, no rigorous arguments to support this belief, and several controversies arose about this subject over the years. It would be evidently challenging to shed more light on this problem using cold atoms in disordered lattices.
3D scaling theory predicts a critical value of disorder, above which every eigenfunction exponentially localizes, and this fact has found strong evidence in numerical simulations.
In the area of AMO physics, effects of disorder have been studied in the context of weak localization of light in random media ͓69͔, which is believed to be a precursor of Anderson localization, and in the form of the so-called dynamical localization, that inhibits diffusion over the energy levels ladder in periodically driven quantum chaotic systems, such as kicked rotor ͓70͔, microwave driven hydrogenlike atom ͑see ͓13͔ and references therein͒, or cold atoms kicked by optical lattices ͓71͔.
It is also worth mentioning at this point the existing large literature on unusual band structure and conductance properties in systems with incommensurate periodic potentials ͓72͔. The famous Harper's equation describing electron's hopping in a cos͑·͒ potential in 1D ͓73͔ may have, depending on the strength of the potential, only localized, or only extended states due to the, so-called, Aubry self-duality. In more complicated cases without self-dual property, and/or in higher dimensions coexistence of localized and extended states is very frequent.

B. Localization in Fermi liquids
The effects of disorder in electronic gases ͑i.e., Fermi gases with repulsive interactions͒ were in the center of interest over many decades. Originally, it was believed that weak disorder should not modify essentially the Fermi liquid quasiparticle picture of Landau. Altshuler and Aronov ͓74͔, and independently Fukuyama ͓75͔, have shown, however, that even weak disorder leads to surprisingly singular corrections to electronic density of states near the Fermi surface, and to transport properties.
As we discussed above, for sufficient disorder in 3D all states are localized, and the standard Fermi liquid theory is not valid. One can use then a Fermi-liquid-like theory using localized quasiparticle states. The system enters then an insulating Fermi glass state ͓76͔, termed often also as Anderson insulator, in which most of the interaction effects are included in the properties of the Landau's quasiparticles.
In 1994 Shepelyansky ͓77͔ stimulated further discussion about the role of interactions by considering two interacting particles ͑TIP͒ in a random potential. He argued that two repulsing or attracting particles can propagate coherently on a distance much larger than the one-particle localization length. Several groups have tried to study these effects of interplay between the disorder and ͑repulsive͒ interactions in more detail in the regime when Fermi liquid becomes unstable as the Mott insulator state is approached by increasing the interactions. Numerical studies performed for spinless fermions with nearest-neighbor ͑NN͒ interactions in a disordered mesoscopic ring; for spin 1 / 2 electrons in a ring, described by the half-filled Hubbard-Anderson model; for spinless fermions with Coulomb repulsion ͑reduced to NN repulsion͒ in 2D, etc. ͓͑78,79͔͒ show that as interactions become comparable with disorder, delocalization takes places. In a 1D ring it leads to the appearance of persistent currents, in 2D the delocalized state exhibits also an ordered flow of persistent currents, which is believed to constitute a novel quantum phase corresponding to the metallic phase observed in experiments, for instance, with a gas of holes in GaAs heterostructures for the similar range of parameters.
Another intensive subject of investigation concerns metal ͑Fermi liquid͒-insulator transition driven by disorder in 3D. Theoretical description of this phenomenon goes back to the seminal works of Efros and Shklovskii ͓80͔ and McMillan ͓81͔. In this context, particularly impressive are the recent results of experiments on disordered alloys, such as amorphous NbSi ͓82͔, where the evidence for scaling and quantum critical behavior was found. Weakly doped semiconduc-tors provide a good model of a disordered solid, and their critical behavior at the metal-insulator transition has been intensively studied ͑cf. ͓83͔͒. Very interesting results concerning in particular various forms of electronic glass: from Fermi glass, with negligible effects of Coulomb repulsion, to Coulomb glass ͓84͔, dominated by the electronic correlations were obtained in the group of Dressel ͓85͔.
Although there exists experimental evidence for delocalization, enhanced persistent currents, and novel metallic phases at the frontier between the Fermi glass and Mott insulator, the further experimental models that physics of cold atoms might provide are highly welcomed. Especially, since the cold atoms Hubbard toolbox should allow one to design with great fidelity the models studied by theorists: spinless fermions extended Hubbard model in 1D, 2D, and 3D, and spin 1 / 2 Hubbard model in a disordered potential, or even more exotic systems such as Fermi systems with SU͑N͒ "flavor" symmetry ͓86͔. Perhaps a less ambitious, but still interesting challenge is to use ultracold atomic gases to create both Fermi glass and a fermionic Mott insulator, and investigate their properties in detail.

C. Localization in Bose systems
At this point it is also worth mentioning the existing literature on the influence of repulsive interactions on Anderson localization in bosonic systems. In the weakly interacting case, one observes at low temperatures the phenomenon of Bose-Einstein condensation ͑BEC͒, but to the eigenstates of the random potential ͑which are Anderson localized͒. Strong nonlinear interactions tend, however, to destroy the localization effects by introducing screening of disorder by the nonlinear mean field potential ͓87,88͔. This happens as soon as the disorder localization length, L dis , becomes larger than the healing length, l heal . Such destruction of localization by weak nonlinearity occurs also in the context of chaos, as discussed in 1993 by Shepelyansky ͓89͔. Several experiments, aiming at observation of localization with BECs have been recently performed with elongated condensates in the presence of a speckle pattern and 1D optical lattices ͓27-29͔. In particular, transport suppression has been observed in the Orsay and Lens experiments, whereas, as we have shown in the Hannover setup ͓29͔, conditions for Anderson localization can be achieved using additional incommensurate superlattices. As the nonlinearity ͑i.e., number of atoms͒ grows the condensate wave function becomes a superposition of exponentially localized modes of comparably low energies. Overlapping of those modes signifies the onset of the screening regime. We believe that similar effects hold in the strongly interacting limit in optical lattices, when they occur at the crossover from the Anderson glass ͑in the weak interacting regime͒ to Bose glass ͑in the strong interacting regime͒ behavior ͑see ͓21͔, and also ͓90͔͒.

D. Localization in superconductors
Obviously, the effects of disorder on superconductivity were studied practically from the very beginning of the theory of superconductors. Already in the late 1950's Anderson and Gorkov considered "dirty" superconductors ͓91͔.
For a weak disorder, Bardeen-Cooper-Schrieffer ͑BCS͒ theory is still valid, but must be modified; the critical temperature is reduced by the localization effects ͓92͔.
The situation is more complex in the case of strong disorder. For example, in 2D superconductors the superconducting state exists only for sufficiently small values of the disorder. This state is often termed a superconducting vortex glass. Cooper pairs in this state condense and form a delocalized "Bose" condensate. This condensate contains a large number of quantum vortices that are immobile and localized in the random potential energy minima associated with disorder ͓93͔. As disorder grows, the system enters the insulating phase, which is a Bose glass of Copper pairs ͑for general theory see ͓26͔͒. Finally, for even stronger values of the disorder the system enters the insulating Fermi glass phase, when the Cooper pairs break down. Obviously, this picture becomes even more complex at the BCS-BEC crossover.
Superconductor-insulator transition has been recently intensively studied in thin metal films on Ge or Sb substrates, that induce disorder on the atomic scale. For not too thin films, transition to superconduting state occurs via Kosterlitz-Thouless-Berezinsky mechanism, whereas for ultrathin films localization effects suppress superconductivity ͓94,95͔. In particular, scaling behavior and scaling exponent were studied in thin bismuth films ͓96͔.
As before, the physics of cold gases might contribute here significantly to our understanding of the influence of quenched disorder on the phenomenon of superconductivity.

E. Localization and percolation
Percolation is a classical phenomenon that is very closely related to localization ͓2͔. Percolation provides a very general paradigm for a lot of physical problems ranging from disordered electric devices ͓97͔, forest fires, and epidemics ͓10,98͔, to ferromagnetic ordering ͓8͔. In lattices, one considers site and bond percolation, and asks the following question: given a probability of filling a lattice site ͑filling a bond͒, and given a layer of the lattice of linear width L, does a percolating cluster of filled sites ͑bonds͒ that connects the walls of the layer exist?
Obviously, a percolating medium with a percolating cluster of empty sites is an example of a medium consisting of randomly distributed scatterers. One has to expect that Anderson localization will occur if quantum waves will propagate and scatter in such a medium. An interplay between percolation and localization has been a subject of intensive studies in recent years. On one hand, when a classical flow is possible, the quantum one might be suppressed due to the destructive interferences and Anderson mechanism. On the other hand, quantum mechanics offers a possibility of tunneling through the classically forbidden regions, and may thus allow for a classically forbidden flow. It turns out that this latter mechanism is very weak, and one typically observes three regimes of localization-delocalization behavior: classical localization below the percolation limit, quantum localization above the percolation limit, and quantum delocalization sufficiently above the percolation limit ͓99,100͔. Quantum percolation plays a role of mechanism responsible for quantum Hall effect ͓101͔. Obviously, interactions in the presence of quantum percolation introduce additional complexity into the phenomenon.
Atomic Fermi-Bose mixtures and atomic gases in general offer an interesting possibility to study quantum percolation in a controlled way. One should stress that quenching atoms as random scatterers in a lattice ͑below percolation thresh-old͒ would be one of the methods itself to generate random local potentials.

F. Random field Ising model
Particularly interesting are those disordered systems, in which arbitrarily small disorder causes large qualitative effects, with Anderson localization in 1D and, most presumably, in 2D being paradigm examples. Other examples are provided by classical systems that exhibit long range order at the lower critical dimension. In such systems, addition of an arbitrary small local potential ͑magnetic field͒, that has a distribution assuming the same symmetry as the considered model, destroys long range order.
The first example of such behavior has been shown by Imry and Ma ͓102͔, using the domain wall argument; it concerns random field Ising model in 2D, for which magnetization vanishes in a random magnetic field in the Ising spin direction with symmetric distribution ͑Z 2 symmetry͒. This result was soon after proven rigorously ͓103͔, and even generalized to XY, Heisenberg, or Potts models ͑provided that the corresponding "field" assumes the same symmetry as the model, i.e., U͑1͒, SU͑2͒, etc. ͓104͔͒.
One should note that most of the above discussed effects concern abstract spin models, and have no direct experimental realizations in condensed matter systems. Cold atoms offer here a unique possibility of both feasible realization of classical models, and of studying quantum effects in those systems. Equally interesting in this context could be spin models in which the random magnetic field breaks the symmetry, such as, for instance, XY model in 2D in the random field directed in, say, the X direction. Such field breaks the U͑1͒ symmetry and changes the universality class of the model to the Ising class. Simultaneously, it prevents spontaneous magnetization in the X direction. In effect, the system attains the macroscopic magnetization in the Y direction ͓105͔. We have recently studied these kinds of systems and proved this result at T = 0 rigorously. We expect in fact finite T transition ͑as in Ising model͒, but a detailed analysis of that case goes beyond the scope of the present paper ͓106͔.

G. Spin glasses: Parisi's theory and the "droplet" model
Spin glasses are spin systems interacting via random couplings, that can be both ferro-or antiferromagnetic. Such variations of the couplings lead typically to frustration. Spin glass models may thus exhibit many local minima of the free energy. For this reason, spin glasses remain one of the challenges of the statistical physics and, in particular, the question of the nature of their ordering is still open. According to Parisi's picture, the spin glass phase consists of very many pure thermodynamic phases. The order parameter of a spin glass becomes thus a function characterizing the probability of overlaps between the distinct pure phases ͓6͔. According to the so-called "droplet" picture, developed by Huse-Fisher and Bray-Moore ͓7,107͔ there are ͑for Ising spin glasses͒ just two pure phases ͑up to Z 2 symmetry͒, and what frustration does is to change very significantly the spectrum of excitations ͑domain walls, droplets͒ close to the equilibrium. While the Parisi's picture ͑related also to the replica symmetry breaking͒ is most presumably valid for long range spin glass models, such as Sherrington-Kirkpatrick model ͓108͔, the "droplet" model is formulated as a scaling theory, and has a lot of numerical support for short range models, such as the Edwards-Anderson ͓109͔ model. For more details of these two pictures relevant for our actual study see Sec. V.

III. DYNAMICS OF COMPOSITE FERMIONS IN THE STRONG COUPLING REGIME
In this section we begin our detailed discussion of the low temperature physics of the Fermi-Bose mixtures. In particular we consider a mixture of ultracold bosons ͑b͒ and spinless ͑or spin-polarized͒ fermions ͑f͒, for example 7 Li-6 Li or 87 Rb-40 K, trapped in an optical lattice. In the following, we will first consider the case of an homogeneous optical lattice, where all lattice sites are equivalent, and we will review previous results focusing on the formation of composite fermions and the quantum phase diagram ͓52͔. Second, we shall extend the analysis to the case of inhomogeneous optical lattices. We consider on-site inhomogeneities consisting of a harmonic confining potential and/or diagonal disorder. In all cases considered below, the temperature is assumed to be low enough and the potential wells deep enough so that only quantum states in the fundamental Bloch band for bosons or fermions are populated. Note that this requires that the filling factor for fermions f , is smaller than 1, i.e., the total number of fermions, N f , is smaller than the total number of lattice sites N.
In the tight-binding regime, it is convenient to project wave functions on the Wannier basis of the fundamental Bloch band, corresponding to wave functions well localized in each lattice site ͓110,111͔. This leads to the Fermi-Bose Hubbard ͑FBH͒ Hamiltonian ͓8,12,19,63͔: and f i are bosonic and fermionic creationannihilation operators of a particle in the ith localized Wannier state of the fundamental band and are the corresponding on-site number operators. The FBH model describes: ͑i͒ nearest-neighbor ͑NN͒ boson ͑fermion͒ hopping, with an associated negative energy, −J b ͑−J f ͒; ͑ii͒ onsite boson-boson interactions with an energy V, which is supposed to be positive ͑i.e., repulsive͒ in the reminder of the paper; ͑iii͒ on-site boson-fermion interactions with an energy U, which is positive ͑negative͒ for repulsive ͑attractive͒ interactions; and ͑iv͒ on-site energy due to interactions with a possibly inhomogeneous potential, with energies − i b and − i f ; eventually, − i b and − i f also contain the chemical potentials in grand canonical description. For the sake of simplicity, we shall focus, in the following, on the case of equal hopping for fermions and bosons, J b = J f = J and we shall assume the strong coupling regime, i.e., V , U ӷ J. Generalization to the case J b J f is just straightforward.

A. Quantum phases in homogeneous optical lattices
Before turning to inhomogeneous optical lattices, let us briefly review here the results presented in ͓52͔ for homogeneous lattices at zero temperature, when all sites are translationally equivalent. In the limit of vanishing hopping ͑J =0͒ with finite repulsive boson-boson interaction V, and in the absence of interactions between bosons and fermions ͑U =0͒, the bosons are in a Mott insulator ͑MI͒ phase with exactly ñ = b + 1 bosons per site, where b = b / V and x denotes the integer part of x. In contrast, the fermions can be in any set of Wannier states, since for vanishing tunneling, the energy is independent of their configuration. The situation changes when the interparticle interactions between bosons and fermions, U, are turned on. In the following, we define ␣ = U / V and we consider the case of a bosonic MI phase with ñ bosons per site. The presence of a fermion in site i may attract −s Ͼ 0 bosons or equivalently expel s ഛ ñ boson͑s͒ depending on the sign of U. The on-site energy gain in attracting −s bosons or expelling s bosons from site i is Within the occupation number basis, excitations correspond to having ñ − s ± 1 bosons in a site with a fermion, instead of ñ − s and, therefore, the corresponding excitation energy is ϳV. In the following, we assume that the temperature is smaller than V so that the population of the abovementioned excitations can be neglected. It follows that tunneling of a fermion is necessarily accompanied by the tunneling of −s bosons ͑if s Ͻ 0͒ or opposed-tunneling of s bosons ͑if s ജ 0͒. The dynamics of our Fermi-Bose mixture can thus be regarded as the one of composite fermions made of one fermion plus −s bosons ͑if s Ͻ 0͒ or one fermion plus s bosonic holes ͑if s ജ 0͒. The annihilation operators of the composite fermions are ͓52͔ These operators are fermionic in the sub-Hilbert space generated by ͉n-ms , m͘ with m =0,1 in each lattice. Note that within the picture of fermionic composites, the vacuum state corresponds to the MI phase with ñ boson per site. At this point, different composite fermions appear depending on the values of ␣, ñ, and b as detailed in Fig. 1 ͓52͔. The different composites are denoted by Roman numbers I, II, III, etc, which denote the total number of particles that form the corresponding composite fermion. Additionally, a bar over a Ro-man number indicates composite fermions formed by one bare fermion and bosonic holes, rather than bosons. For the sake of simplicity, we shall consider only bosonic MI phases with ñ = 1 boson per site ͑i.e., 0 Ͻ b ഛ 1͒ in the following parts of this paper ͓113͔. If ␣ − b Ͼ 0, a fermion in site i pushes the boson out of the site. We will call the sites with this property B-sites. This notation will become particularly important in the presence of disorder ͑local b ͒. The composites, in this case, correspond to one fermion plus one bosonic hole ͓this phase is called II in Fig. 1͑a͔͒. If −1Ͻ ␣ − b Ͻ 0, we have bare fermions as composites ͑this corresponds to phase I͒. The sites with this property will be called A-sites. Finally, if −2 Ͻ ␣ − b Ͻ −1, the composites are made of one fermion plus one boson ͑phase II͒. The sites with the latter property will be called C-sites. Because all sites are equivalent for the fermions, the ground state is highly ͓N!/͑N f ͒!͑N − N f ͒!͔-degenerated, so the manifold of ground states is strongly coupled by fermion or boson tunneling. We assume now that the tunneling rate J is small but finite. Using time-dependent degenerate perturbation theory ͓114͔, we derive an effective Hamiltonian ͓12͔ for the fermionic composites:

͑4͒
where M i = F i † F i and eff is the chemical potential, whose value is fixed by the total number of composite fermions. The nearest neighbor hopping for the composites is described by −d eff and the nearest neighbor compositecomposite interactions are given by K eff , which may be repulsive ͑Ͼ0͒ or attractive ͑Ͻ0͒. This effective model is equivalent to that of spinless interacting fermions. The interaction coefficient K eff originates from second order terms in perturbative theory and can be written in the general form:

͑5͒
This expression is valid in all the cases but when s = 0 the last term ͑1/s␣͒ should not be taken into account. d eff originates from ͉͑s͉ +1͒-th order terms in perturbative theory and thus presents different forms in different regions of the phase diagram of Fig. 1. For instance in region I, d eff = J, in region II d eff =2J 2 / U, and in region II, d eff =4J 2 / ͉U͉. The physics of the system is determined by the ratio K eff / d eff and the sign of K eff . In Fig. 1͑a͒, the subindex A/R denotes attractive ͑K eff Ͻ 0͒/repulsive ͑K eff Ͼ 0͒ composites interactions. Figure 1͑b͒ shows the quantum phase diagram of composites for fermionic filling factor f = 0.4 and tunneling J / V = 0.02. As an example, let us consider the case of repulsive interactions between bosons and fermions, ␣ Ͼ 0. Once the fermion-bosonic hole composites II have been created ͑␣ Ͼ b ͒, the relation K eff / d eff =−2͑␣ −1͒ applies. Consequently, if one increases the interactions between bosons and fermions adiabatically, the system evolves through different quantum phases. For b Ͻ ␣ Ͻ 1, the interactions between composite fermions are repulsive and of the same order of the tunneling; the system exhibits delocalized metallic phases. For ␣ Ӎ 1, the interactions between composite fermions vanish and the system shows up the properties of an ideal Fermi gas. Growing further the repulsive interactions between bosons and fermions, the interactions between composite fermions become attractive. For 1 Ͻ ␣ Ͻ 2, one expects the system to show superfluidity, and for ␣ Ͼ 2 fermionic insulator domains are predicted.
In the reminder of the paper, we shall generalize these results to the case of inhomogeneous optical lattices. We shall assume diagonal inhomogeneities, i.e., site-dependent local energies ͑ i b,f depends on site i but the tunneling rate J and interactions U and V do not͒. Diagonal inhomogeneities may account for ͑i͒ overall trapping potential ͑usually har-monic͒, which is usually underlying in experiments on ultracold atoms, and ͑ii͒ disorder that may be introduced in different ways in ultracold samples ͑see Sec. VI for details͒.

B. Composites in disordered lattices: Effective Hamiltonian
In this section we include on-site energy inhomogeneities in the optical lattice and we derive a generalized effective Hamiltonian for the composite fermions. Strictly speaking, in the presence of disorder the hopping terms should depend on the site under consideration. Nevertheless, for weak enough disorder one can assume site independent tunneling for both bosons or fermions ͓21͔. In the following we will restrict ourselves to the case where the hopping rates of bosons and fermions are equal and site-independent, J b = J f = J and to the strong coupling regime, V , U ӷ J, where the tunneling can be considered as a perturbation, as in Sec. III A.
For homogeneous lattices ͑see Sec. III A͒, following the lines of Refs. ͓51,52͔, we have used the method of degenerate second order perturbation theory to derive the effective Hamiltonian ͑4͒ by projecting the wave function onto the multiply degenerated ground state of the system in the absence of tunneling.
In the inhomogeneous case, this approach cannot be applied since even for J = 0 there exists a well-defined single ground state determined by the values of the local chemical potentials. Nevertheless, in general, there will be a manifold of many states with similar energies. The differences of energy inside a manifold are of the order of the difference of chemical potential in different sites, whose random distribution is bounded, i.e., 0 ഛ i b , i f ഛ 1. Moreover, the lower energy manifold is separated from the exited states by a gap given by the boson-boson interaction, V. Therefore one can apply a form of quasidegenerate perturbation theory by projecting onto the manifold of near-ground states ͓114͔.
As it is described in Ref. ͓114͔ and briefly summarized in Appendix A, we construct an effective Hamiltonian that describes the slow, low-energy perturbation induced within the manifold of unperturbed ground states by means of a unitary transformation applied to the total Hamiltonian ͑1͒. By denoting with P the projector on the manifold and Q =1− P its complement, the expression of the effective Hamiltonian can be written as As second order theory can only connect states that differ on one set of two adjacent sites, the effective Hamiltonian H eff can only contain nearest-neighbor hopping and interactions as well as on-site energies i ͓37͔: where M i , F i are defined as in Eq. ͑4͒. The explicit calculation of the coefficients d ij , K ij , and i depends on the concrete type of composites. In the three following sections we address the cases of fermion-bosonic hole composites ͑II͒, bare fermion composites ͑I͒, and fermion-boson composites ͑II͒.

C. Fermion-Bosonic hole composites
In this section, we assume that all sites are B-sites, i.e., ␣ − i b Ͼ 0, so that composite fermions II are created. This means that each site contains either one boson or one fermion plus a bosonic hole. Thus the manifold of low lying states comprises all possible configurations of N f fermions on an N-site lattice, with no fermion occupied sites filled by bosons.
Within the manifold of ground states, a fermion jump from site i to site j can only occur if the boson that was initially in site j jumps back to site i into the hole the fermion leaves behind. Therefore the number operator for fermions and bosons are related with the number operator of composites, i.e., M i = m i =1−n i . Note that the composite model is expressed in terms of the composite fermionic operators To determine the coefficients in Eq. ͑7͒, one looks at two adjacent sites with indices i and j and uses a vector notation ͉1 b ,1 f ͘ which would correspond to one boson on site i and one fermion on site j. In the composite fermion picture this would be denoted as ͉0,1͘ c , i.e., one composite fermion on site j and no composite fermion on site i. With this notation tunneling rates and nearestneighbor interactions are calculated from Eq. ͑6͒ as Summing these terms up yields the coefficients for Eq. ͑7͒: Here, ͗i , j͘ represents all neighbor sites of i. We shall now consider separately two limiting cases: ͑i͒ i f = 0 and i b = i V, and ͑ii͒ i f = i b = i V.

Case where i f =0
In the first case, we assume that the on-site energy for fermions vanishes. We assume also that all sites are B-sites, i.e., ␣ − i b Ͼ 0 everywhere. In this case, the hopping amplitudes d ij are always positive, although may vary quite significantly with disorder, especially when ⌬ ij b Ӎ ␣. As shown in Fig. 2, for ␣ Ͼ 1, K ij ഛ 0 and we deal with attractive ͑although random͒ interactions. For ␣ Ͻ 1, K ij ജ 0 and the interactions between composites are repulsive. For ␣ Ͻ 1, but close to 1, K ij might take positive or negative values for ⌬ ij b small or ⌬ ij b Ӎ ␣. In this case, the qualitative character of interactions may be controlled by inhomogeneity ͓37͔. At low temperatures the physics of the system will depend on the relation between i b 's and ␣. (a) Small disorder limit. For small disorder, we may neglect the contributions of ⌬ ij b to d ij Ӎ d and K ij Ӎ K, and keep only the leading disorder contribution in i , i.e., the first term in Eq. ͑14͒. Note that the latter contribution is relevant in 1D and 2D leading to Anderson localization of single particles ͓66͔. When K / d Ӷ 1 the system will then be in the Fermi glass phase, i.e., Anderson localized ͑and many-body corrected͒ single particle states will be occupied according to the Fermi-Dirac rules ͓76͔. For repulsive interactions and K / d ӷ 1, the ground state will be a Mott insulator and the composite fermions will be pinned for large filling factors. In particular, for filling factor f =1/2, one expects the ground state to be in the form of a checkerboard. For intermediate values of K / d, with K Ͼ 0, delocalized metallic phases with enhanced persistent currents are possible ͓78͔. Similarly, for attractive interactions ͑K Ͻ 0͒ and ͉K͉ / d Ͻ 1 one expects competition between pairing of fermions and disorder, i.e., a "dirty" superfluid phase while for ͉K͉ / d ӷ 1, the fermions will form a domain insulator. Figure 3 shows a schematic representation of expected disordered phases of the type II fermionic composites for small disorder, and vanishing fermionic on-site chemical potential.
(b) Spin glass limit. Another interesting limit corresponds to the case ⌬ ij b Ӎ ␣ Ӎ 1. Such a situation can be achieved by combining a superlattice potential with a spatial period twice as large as the one of the lattice ͑which alone induces ͉⌬ ij b ͉ =1͒ and a random potential to induce site-to-site fluctuations. The tunneling becomes then nonresonant and can be neglected in Eq. ͑7͒, while the couplings K ij fluctuate strongly as shown in Fig. 2. We end up then with the ͑fermi-onic͒ Ising spin glass model ͓37͔ described by the Edwards-Anderson model with s i =2M i − 1 = ± 1. This case is studied in more detail in Sec. V.

Case where
Let us now consider that the chemical potential is equal for bosons and fermions at each lattice site, i f = i b = i V. All sites are still assumed to be B-sites.
The effective interactions are for ␣ Ͼ 1 always negative, and therefore the composites experience random attractive interactions ͑as in the previous case͒, while for ␣ Ͻ 1, K ij Ͼ 0, and therefore we deal with random repulsive interactions. For ␣ = 1, the interactions between composites vanish for all the values of the amplitude of the disorder.
In this case the sign of the interactions between composites is governed by the interactions between bosons and fermions alone. Note that this is not possible when one considers only disorder for the bosons. Figure 4 shows the tunneling and the nearest-neighbor couplings for different values of ␣. We expect here the appearance of similar phases, as in the previously discussed case.

D. Bare Fermion composites
In this section we now assume that all sites are A-sites and correspond to type I fermionic composites, i.e., −1 Ͻ ␣ − i b Ͻ 0. This means that composite fermions reduce to bare fermions ͑F i = f i ͒ flowing on the top of a MI phase with ñ =1 boson per site. Each site now contains one boson plus eventually one fermion. From application of perturbation theory as described in Sec. III B ͓see Eq. ͑6͔͒, one finds that the

͑17͒
We observe that the inhomogeneities for fermions ͑sitedependent i f ͒ neither perturb the effective tunneling, nor the effective interaction parameter, while i Ӎ − i f up to corrections of the order of O͑J 2 / V͒ for type I composite ͑bare͒ fermions. In this case, composite tunneling d ij originates from the first order term, while the nearest-neighbor interaction originates from second order perturbation. It should be noted that in the case of type I composites, the hopping d ij and interaction K ij parameters in Eq. ͑7͒ do not depend on the sign of the fermion-boson interaction ␣.
The couplings K ij are always positive, and for ␣ Ӎ 0, K ij Ӎ O͑␣ 2 ͒, and both the repulsive interactions and disorder are very weak, leading to an almost ideal Fermi liquid behavior at low temperature. For finite ␣, and ⌬ ij b Ӎ 1−͉␣͉, however, the fluctuations of K ij might be quite large as shown in Fig. 5. Note that for ͉␣͉Ӎ1, this will occur even for small disorder. It is interesting to note that the dynamics of type I composites in our system resembles quantum bond percolation. As suggested from Fig. 5, one can assume in a somehow simplified view that the interaction parameter K ij takes either very large, or zero values. The lattice decomposes into two sublattices ͑see Fig. 6͒: a "weak" bond sub-lattice ͑corresponding to K ij Ӷ d ij ͒ in which fermions flow as in an almost ideal Fermi liquid, and a "strong" sublattice ͑corresponding to K ij ӷ d ij ͒, where only one fermion per bond is allowed ͑M i M j = 0 for all nearest neighbors in the "strong" cluster͒. Therefore we see that the physics of bond percolation ͓10,97͔ will play a role. For p Ͼ p c , where p is the density of weak bonds and p c Ӎ 0.50 ͑in two-dimensional square lattices͒ and p c Ӎ 0.25 ͑in three-dimensional cubic lat-tices͒, the weak bond sublattice will be percolating, i.e., there exists a large cluster of weak bonds which spans the lattice from one side to the other. The question arises as to determine the quantum bond percolation threshold p Q , i.e., for which minimal value of p, the eigenstates of the quantum gas will be delocalized over the extension of the system. Although it is clear that p Q Ͼ p c , it is still an open question to determine the exact value for the quantum percolation threshold p Q ͓99,115-118͔. Therefore experimental realization of our system may be of considerable interest for addressing this general question.

E. Fermion-Boson composites
We finally consider in this section the case, when all sites are C-sites, so that type II composites corresponding to −2 Ͻ ␣ − i b Ͻ −1 are formed. The composites are made of one fermion and one boson. This means that each lattice site is populated by either one boson or one fermion plus two bosons. Tunneling as well as nearest-neighbor interaction of composites arise from second order terms in perturbative theory ͓see Sec. III B and Eq. ͑6͒ for details͔. Along the lines of Sec. III B, we find the following expressions:

͑20͒
Different scenarios also arise in this case. In the following, we shall consider the case i f = 0 and i b = i V. The other extreme case, i f = i b = i V, leads to qualitatively similar conclusions.

Case where i f =0
We assume here that the on-site energy for fermions is i f = 0. As for Fig. 2, we plot the effective tunneling and interaction parameter versus inhomogeneity parameter ⌬ ij b in Fig. 7. The general behavior of d ij and K ij is qualitatively the same as in the case of type II composites. For type II composites, and for small disorder, we find K / d =1−4͉␣͉ +3/͑2 − ͉␣͉͒ with 1 Ͻ ͉␣͉ Ͻ 2. The inhomogeneity is now given, for type II composites, by i =− i b . The regimes, where K Ӷ d corresponding to an almost ideal Fermi gas ͑in the absence of disorder͒, or to a Fermi glass ͑in the presence of disorder͒, can be reached in the region ␣ Ӎ −5 / 4. The opposite regimes of strong effective interactions, where K ӷ d appears for ␣ տ −2 and corresponds to repulsive interactions K Ͼ 0. In this region, the fermionic checkerboard phase if the filling factor is 1 / 2 ͑for vanishing disorder͒ and the repulsive Fermi glass phase ͑in the presence of disorder͒ are expected. Here, no strong attractive interactions regime occurs since K / d reaches a minimum of Ӎ−0.07 for ␣ = ͱ 3/2−2. Therefore, in contrast to type II composites, for type II composites: ͑i͒ due to weakness of attractive interactions, the domain insulator phase does not appear, and even the "dirty" superfluid phase may be washed out; and ͑ii͒ arbitrary strong repulsive interactions can be used to generate a Mott insulator, which might be difficult for II composites, where K / d is limited to 2, suggesting that Fermi liquid, Fermi glass behavior will prevail.
As shown in Fig. 7, the spin glass limit can also be reached, for example, for ␣ Ӎ −1.05 and ͉⌬ ij f ͉Ӎ0.9. In this regime, the tunneling is nonresonant due to strong disorder and the nearest-neighbor interaction fluctuates strongly from negative to positive values. See Sec. V for further study of the spin glass limit.

Sites A and B
Obviously, the situation becomes much more complex when we deal with different types of sites in the lattice. Again there are infinitely many possibilities, and the simplest ones are, for instance: ͑i͒ coexistence of A-and B-sites, or ͑ii͒ A-and C-sites, or ͑iii͒ A-, B-, and C-sites, etc. In the following we shall consider only the case ͑i͒ with i f = 0 and i b = i V, since the other cases lead to qualitatively similar effects.
Let us assume that the numbers of A-and B-sites are macroscopic, i.e., of the order of N. More precisely, we will consider that N A ͑number of A-sites͒ and N B ͑number of B-sites͒ of order N / 2. In this case the physics of site percolation ͓10͔ will play a role. If N f ഛ N B the composite fermions will move within a cluster of B-sites. When N B will be above the classical percolation threshold, this cluster will be percolating. The expressions Eq. ͑12͒ and Eq. ͑13͒ will still be valid, except that they will connect only the B-sites. The physics of the system will be similar as in the case of type I composites͒, but it will occur now on the percolating cluster. For small disorder, and K / d Ӷ 1, the system will be in a Fermi glass phase in which the interplay between the Anderson localization of single particles due to fluctuations of i b and quantum percolation effects, that is randomness of the B-sites cluster, will occur. For repulsive interactions and K / d ӷ 1, the ground state will be a Mott insulator on the cluster and the composite fermions will be pinned ͑in particular for half-filling of the cluster͒. It is an open question whether the delocalized metallic phases with enhanced persistent current of the kind discussed in Ref. attractive interactions ͑K Ͻ 0͒ and ͉K͉ / d Ͻ 1 pairing of ͑perhaps localized͒ fermions will take place. In the case of ͉K͉ / d ӷ 1, we expect that the fermions will form a domain insulator on the cluster. In the spin glass limit ⌬ ij b Ӎ ␣ Ӎ 1, we will deal with the Edward-Anderson spin glass on the cluster. Such systems are of interest in condensed matter physics ͑cf. ͓11͔͒, and again questions connected to the nature of spin glass ordering may be studied in this case.
When N f Ͼ N B , all B-sites will be filled, and the physics will occur on the cluster of A-sites. For ␣ Ӎ 0, we will deal with a gas with very weak repulsive interactions, and no significant disorder on the random cluster; this is an ideal test to study quantum percolation at low T. For finite ␣, and ⌬ ij b Ӎ 1−␣, the interplay between the fluctuating repulsive K ij 's and quantum percolation might be studied.

A. Numerical method
In this section, we present numerical results that give evidences of ͑i͒ formation of composite particles in Fermi-Bose mixtures in optical lattices and ͑ii͒ existence of different quantum phases for various sets of composite tunneling and interaction parameters and inhomogeneities. We mainly focus on type II composites. Mean-field theory provides appropriate although not exact properties of Hubbard models ͓8͔. In the following, we consider a variational mean field approach provided by the Gutzwiller ansatz ͑GA͒ ͓26,119͔. In particular, the GA ansatz has been successfully employed for bosonic systems to study the superfluid to Mott insulator transition ͓19,26͔ in nondisordered lattices, and the Anderson and Bose glass transitions in the presence of disorder ͓21,26͔.
Briefly, the Gutzwiller approach neglects site-to-site quantum coherences so that the many-body ground state is written as a product of N states, each one being localized in a different lattice site. Each localized state is a superposition of different Fock states ͉n , m͘ i with exactly n bosons and m fermions on the ith lattice site: where n max is an arbitrary maximum occupation number of bosons in each lattice site ͓120͔.
The g n,m ͑i͒ are complex coefficients proportional to the amplitude of finding n bosons and m fermions in the ith lattice site, and consequently we can impose, without loss of generality, these coefficients to satisfy ͚ n ͚ m ͉g n,m ͑i͒ ͉ 2 = 1. For the sake of simplicity, we neglect the anticommutation relation of fermionic creation ͑f i ͒ and annihilation ͑f j † ͒ operators in different lattice sites. However, Pauli principle applies in each lattice site ͑m i ഛ 1 ∀ i͒. Since GA neglects correlations between different sites, this procedure is expected to be safe and is commonly used within the Gutzwiller approach ͓53͔.
Inserting ͉ MF ͘ in the Schrödinger equation with the twospecies Fermi-Bose Hamiltonian ͑1͒, we were able to deter-mine the ground state and to compute the dynamical evolution of the Fermi-Bose mixture.

Ground state calculations
Employing a standard conjugate-gradient downhill method ͓121͔, we minimize the total energy ͗ MF ͉H FBH ͉ MF ͘ with H FBH given by Eq. ͑1͒ under the constraint of fixed total numbers of fermions N f and bosons N b ͓122͔: The numerical procedure is as follows: ͑i͒ We minimize the energy of the mixture ͑eventually͒ in the presence of smooth trapping potentials and with nonzero tunneling for bosons and fermions, but assuming vanishing interactions between bosons and fermions ͑U =0͒. During the minimization the normalization ͚͑ n ͚ m ͉g n,m ͉ 2 =1͒ should be imposed. ͑ii͒ After this initial minimization, we ramp up adiabatically the interactions between bosons and fermions using the dynamical Gutzwiller approach ͑see below͒. In this way, we end up with the ground state of the mixture in the presence of tunneling J, nonvanishing interactions U and V, and eventually in the presence of a smooth trapping potential.
This two-step procedure is indeed necessary because in the presence of interactions between bosons and fermions finite numbers of bosons and fermions correspond to a saddle point of Eq. ͑22͒, and no true minimum can be found within direct minimization of the total Hamiltonian ͓54͔.

Time-dependent calculations
Using the time-dependent variational principle ͑͗ MF ͉iប‫ץ‬ t − H FBH ͑t͉͒ MF ͘ → min͒ with Hamiltonian H FBH given by Eq. ͑1͒ and eventually time-dependent parameters J f,b , U, V, i f,b , we end up with the following dynamical equation for the Gutzwiller coefficients ͓54,123͔: ͑j͒ ͬ .

͑25͒
Note that these equations are valid under the hypothesis of neglecting anticommutation relations for fermionic operators in different sites. Equations ͑23͒-͑25͒ preserve both normalization of the wave function and the mean particle numbers. In the following, the dynamical Gutzwiller approach will be used for ͑i͒ computing the ground state of the mixture in the presence of interactions between bosons and fermions ͑see above͒ and ͑ii͒ to ramp up adiabatically disorder in the optical lattice potential.

B. Numerical results
We have considered a 2D optical lattice with N =10ϫ 10 sites to perform the simulations of the different quantum phases that appear for type IĪ composites in the presence of a very shallow harmonic trapping potential ͓ i b,f = b,f ϫ l͑i͒ 2 , where l͑i͒ is the distance from site i to the central site in cell size units͔, with different amplitudes for bosons and fermions. The harmonic on-site energy simulates shallow magnetic or optical trapping. The confining potential is experimentally of vital importance in order to see Mott insulator phases that require commensurate filling ͓20,124,125͔. It plays the role of a local chemical potential, and it has been predicted that it modifies some properties of strongly correlated phases ͓126͔. Additionally, this breaks the equivalence of all lattice sites and makes it more obvious the different phases that one can achieve ͑see below͒. We first calculate the ground state of the system considering N b = 60 bosons, N f = 40 fermions, J b / V = J f / V = 0.02, U = 0 in the presence of harmonic traps characterized by b =10 −7 and f =5ϫ 10 −7 .
Under these conditions, we find that, as expected, the bosons are well inside the MI phase with ñ = 1 boson per site ͓19,26͔. Due to the very small values of f and b neither the bosons nor the fermions feel significantly the confining trap as shown in Fig. 8͑a͒.

Nondisordered phases
Starting with this ground state we adiabatically grow the repulsive interactions between bosons and fermions, U, keeping the repulsion between bosons, V, constant, i.e., growing effectively ␣ in order to create the composites. Once the composites appear, the only nonzero probabilities are ͑i͒ ͉g n=1,m=0 ͑i͒ ͉ 2 to have one boson and zero fermion ͑i.e., no composite͒, or ͑ii͒ ͉g n=0,m=1 ͑i͒ ͉ 2 to have zero boson and one fermion ͑i.e., one composite͒. This proves the formation of type IĪ composites ͓128͔. In Fig. 8͑b͒, we show the probability of having one composite ͉g n=0,m=1 ͑i͒ ͉ 2 ͑i.e., one fermion and zero boson͒ at each lattice site for ␣ = 0.5, which corresponds to repulsive interactions between composites K eff = d eff = 1.6 ϫ 10 −3 . Due to the important value of the composite tunneling d eff , the ground state is delocalized and corresponds to a ͑nonideal͒ Fermi liquid.
Increasing further the fermion-boson interaction parameter, ␣, the system reaches the point where the interactions between composites are negligible corresponding to the region of an ideal Fermi gas phase ͑␣ Ӎ 1͒. Figure 8͑c͒ displays the probability of having a composite in each lattice site in the case where the interactions between composites exactly vanish, i.e., ␣ = 1. Increasing again the interaction parameter ␣, one reaches for ␣ Ͼ 1 the region where the interactions between composites are attractive ͑K eff Ͻ 0͒. In this region, composite fermionic insulator domains are predicted. Due to the attractive interactions, the probability of having composite fermions in the center of the trap increases reaching nearly one for high enough effective attractive interactions as shown in Fig. 8͑d͒.
It is also worth noticing that the energies involving the composites are at least three orders of magnitude smaller than the corresponding energies for bosons and fermions ͑d eff , K eff Ӷ J , U , V͒. As a consequence, the effect of inhomogeneities is much larger for composites than it is for bare bosons and fermions. This is exemplified in Fig. 8. For no interaction between bosons and fermions ͑␣ =0͒, the bare particles are not significantly affected by the harmonic trap on the 10ϫ 10 lattice ͓see Fig. 8͑a͔͒. On the contrary, as soon as composites IĪ are created, the harmonic trapping clearly reflects in inhomogeneous population of the lattice sites ͓see Figs. 8͑b͒-8͑d͔͒. Another important consequence is that large time scales are necessary in time evolution processes in order to fulfill the adiabaticity condition.

Disordered phases
We now consider disordered optical lattices for the bosons. The on-site energy i b is assumed to be random with and independent from site to site. For this, we create type IĪ composites in different regimes ͑this is controlled by the value of ␣ as shown before͒ and we slowly ramp up the disorder from 0 to its final value ⌬. Let us first consider a Fermi gas in the absence of disorder ͓see fermions become more and more localized in the lattice sites to form a Fermi glass. Indeed, Fig. 9͑a͒ shows that the fluctuations in composite number are significantly reduced as the amplitude of the disorder increases. For ⌬ =5ϫ 10 −4 , the composite fermions are pinned in random sites as shown in Fig. 9͑c͒. As expected, the N f composite fermions populate the N f sites with minimal i b . It should be noted that in the absence of interactions between bosons and fermions ͑i.e., when the composites are not formed͒, no effect of disorder is observed. This again shows the formation of composites with typical energies significantly smaller than those of bare particles.
We now consider the Fermi insulator domain phase ͓see Fig. 10͑b͔͒ with slowly increasing disorder. In the Fermi insulator domain ͑in the absence of disorder͒, the composite fermions are pinned in the central part of the confining potential. In addition, there is a ring of delocalized fermions and this gives finite fluctuations on the atom number per site ͓ͱ͗͑m i − ͗m i ͒͘ 2 ͘Ӎ0.35͔. As shown in Fig. 10͑a͒, while ramping up the amplitude of disorder, the fluctuations decrease fast and reach ͱ͗͑m i − ͗m i ͒͘ 2 ͘Ӎ0.13 for ⌬Ͼ10 −4 . This indicates that the composite fermions are pinned in different lattice sites. This is confirmed in Fig. 10͑c͒ where we plot the population of the composites fermions in each lattice site for ⌬ =5ϫ 10 −4 . Contrary to what happens for the transition from Fermi gas to Fermi glass, the composites mostly populate the central part of the confining potential. The reason for that is twofold. First, with our parameters, the attractive interaction between composites is of the order of K Ӎ −1.4 ϫ 10 −3 and can compete with disorder ⌬ =3ϫ 10 −4 . This explains the central insulator domain. Second, because tunneling is small ͑d Ӎ 8 ϫ 10 −5 ͒ and because disorder breaks the symmetry of lattice sites in the ring around the domain, the atoms in this region get pinned. The populated sites match the lowest i b .

V. SPIN GLASSES
In this section we discuss in more detail the possible realization of the Edwards-Anderson spin glass Fermi-Bose mixtures as discussed in Sec. III C. Strictly speaking, since the system is quantum it allows for realization of fermionic spin glass ͓129͔. The main goal of such an investigation is to study the nature of the spin glass ordering and to compare the predictions of the Mézard-Parisi and "droplet" pictures.
Although we work along the lines of the original papers ͓6͔, it is necessary to reformulate the standard Mézard-Parisi mean field description of our system. The main difference appears because the Ising spins are coded as presence or absence of a composite at a given site. This leads to a fixed magnetization due to the fixed number of particles in the system. For this we repeat very shortly the Sherrington-Kirkpatrick calculations ͓108͔ here, adapted to our case. FIG. 9. Dynamical crossover from the Fermi gas to the Fermi glass phases. The parameters are the same as in Fig. 8͑c͒. ͑a͒ Variance of the number of fermions per lattice site as a function of the amplitude of the disorder ⌬. ͑b͒ Probability of having one composite ͑one fermion and zero boson͒ at each lattice site for the N sites in the absence of disorder and ͑c͒ after ramping up adiabatically diagonal disorder with amplitude ⌬ =5ϫ 10 −4 .
FIG. 10. Dynamical crossover from the fermionic domain insulator to a disordered insulating phase. The parameters correspond to Fig. 8͑d͒. ͑a͒ Variance of the number of fermions per lattice site as a function of the amplitude of the disorder. ͑b͒ Probability of having one composite ͑one fermion and zero boson͒ in each lattice site in the absence of disorder and ͑c͒ after ramping up adiabatically the disorder with amplitude ⌬ =3ϫ 10 −4 .

A. Edwards-Anderson model for composite fermions
The spin glass limit obtained in Sec. III C with large disorder derives of the composite fermionic model ͑7͒ with vanishing hopping due to strong site-to-site energy fluctuations and NN interactions, K ij . By appropriate choice of ⌬ ij b , K ij fluctuate around mean zero with random positive and negative values ͓see Fig. 2͑b͔͒. Replacing the composite number operators with a classical Ising spin variable s i ª 2M i − 1 = ± 1, one ends up with the Hamiltonian: It describes an ͑fermionic͒ Ising spin glass ͓129͔, which differs from the Edwards-Anderson model ͓6,130͔ in that it has an additional random magnetic field i and, moreover, has to satisfy the constraint of fixed magnetization value, m =2N f / N − 1, as the number of fermions in the underlying FBH model is conserved. It, however, shares the basic characteristics with the Edwards-Anderson model as being a spin Hamiltonian with random spin exchange terms K ij . In particular, this provides bond frustration, which in this model is essential for the appearance of a spin glassy phase. The experimental study of this limit thus could present a way to address various open questions of spin glass physics concerning the nature and the ordering of its ground and possibly metastable states ͑the Mézard-Parisi picture ͓6͔ versus the "droplet" picture ͓7,107͔͒, broken symmetry, and dynamics in classical ͑in the absence of hopping͒ and quantum ͑with small, but nevertheless present hopping͒ spin glasses ͓8,131͔.
For sufficiently large systems, Eq. ͑26͒ is well approximated by assuming K ij and i to be independent random variables with Gaussian distribution, with mean 0 and H, respectively, and variances K / ͱ N and h, respectively ͓132͔. This approximation will be used in the following calculations.
Before employing the mean field approach for Edwards-Anderson-like models in Sec. V E for the Hamiltonian ͑26͒, a very basic outline of the different phases of the short-range Ising models with bond frustration is given and the two competing physical pictures for the spin glass phase are briefly summarized in this section.
The experimental observations have led to the identification of three equilibrium phases, which are characterized by two order parameters ͑for zero external magnetic field͒: M ª ͗s i ͘ T is the magnetization, i.e., the order parameter for magnetic ordering, and Q EA ª ͗s i ͘ T 2 is the Edwards-Anderson order parameter for spin glass ordering. Here, ͗·͘ T denotes the Gibbs ensemble average and · the disorder average. The three phases are ͑i͒ an unordered paramagnetic phase, with M = 0 and Q EA =0; ͑ii͒ an ordered spin glass phase with M = 0 and Q EA 0 that is separated from the paramagnetic phase by a second order phase transition ͓133͔; and ͑iii͒ dependent on the mean value of K ij , an ordered ferromagnetic phase with M 0 and Q EA = 0. It should be pointed out that there are additional questions-different from, but of course connected to the ones discussed in the following-about the nature of the equilibrium spin glass state that stem from the intrinsic problems that are associated in this system with separating equilibrium from nonequilibrium effects such as metastability, hysteresis, and others ͑see ͓134͔ and references therein͒.

B. Mézard-Parisi picture
The Mézard-Parisi ͑MP͒ picture is fundamentally guided by the results of the mean-field theory. On the level of statistical physics, the Gibbs equilibrium distribution of the spin system in the MP picture at temperature T and for a particular disorder configuration K can be written as a unique convex combination of infinitely many pure equilibrium state distributions ͓7,135͔, where the overlap between two pure states is defined as where ⍀ denotes the size of the system. The mean-field version of Q ␣␤ emerges naturally from the calculation in the next section and motivates definition ͑28͒.
In the MP picture, the spin glass transition is interpreted akin to the transition from an Ising paramagnet to a ferromagnet. There, the Gibbs distribution is written as a sum of only two pure states, corresponding to the two possible fully spin-polarized ferromagnetic ground states. As the temperature of the system decreases, the Z 2 symmetry of the system is broken, and a phase transition to a ferromagnetic phase occurs, whose equilibrium properties are not described by the Gibbs state, but by the relevant pure state distribution alone ͓6͔. Analogously, the spin glass transition is characterized by the breaking of the infinite index symmetry, called replica symmetry breaking in the mean-field case, by which one pure state distribution T,K ␣ is chosen and alone describes the low-temperature properties of the system ͓6͔. However, unlike the Ising ferromagnet, the pure states of the spin glass are not related to each other by a symmetry of the Hamiltonian, but rather by an accidental, infinite degeneracy of the ground state caused by the randomness of the bonds and the frustration effects. This picture can be interpreted as the system getting frozen into one particular state out of infinitely many different ground or metastable states of the system. These states are all taken to be separated by free energy barriers, whose height either diverges with the system size or it is finite but still so large that the decay into a "true" ground state does not occur on observable time scales. Thus fluctuations around one of these ground states can only sample excited states within one particular free energy valley. Consequently, Q EA in the spin glass phase must be redefined as Q EA ª Q ␣␣ , the self-overlap of the state, whereas it remains unchanged for the paramagnetic phase.

C. de Almeida-Thouless plane
Based on the results of mean-field theory, one of the predictions of the MP picture concerns the order of the infinitely many spin glass ground states, which is ultrametric ͓6,136͔, as can be seen from the joint probability distribution of three different ground state overlaps, P K ͑Q 12 , Q 13 , Q 23 ͒. Upon choosing independently three pure states 1, 2, and 3 from the decomposition ͑27͒, one should find that with probability 1/4, Q 12 = Q 13 = Q 23 and with probability 3 / 4 two of the overlaps are equal and smaller than the third. Ultrametricity then follows from the canonical distance function D ␣␤ = Q EA − Q ␣␤ . The mean-field theory, both with and without a magnetic field, predicts the existence of a plane in the space of the Hamiltonian parameters, called the de Almeida-Thouless ͑dAT͒ plane ͓137͔, below which the naive ansatz for the spin glass phase becomes invalid and the system is characterized by the transition to this ultrametrically ordered infinite manifold of ground states. It should be pointed out that the clear occurrence of such a dAT plane in the finite range Edwards-Anderson model would be an important indicator for the validity of the MP picture in these systems. As we discuss in Sec. V D, this conclusion has to be drawn with great care.

D. "Droplet" model
The very applicability of the MP model for finite-range systems is, however, still unproven. It is both challenged by a rival theory, the so-called droplet model ͓107͔, as well as by mathematical analysis ͑cf. ͓7͔ and references therein͒ that questions the validity of transferring a picture developed for the infinite-range mean-field case to the short-range model. Being a phenomenological theory based on scaling arguments and numerical results, the droplet model describes the ordered spin glass phase below the transition as one of just two possible pure states, connected by spin-flip symmetry, analogous to the ferromagnet mentioned above. Consequently there can be no infinite hierarchy of any kind, and thereby no ultrametricity. Excitations over the ground state are regions with a fractal boundary-the droplets-in which the spins are in the configuration of the opposite groundstate. The free energy of droplets of diameter L is taken to scale as ϳL , with Ͻ 0 at and below the critical dimension, which is generally taken to be two. So three is the only physical dimension where the spin glass transition is stable with nonzero transition temperature, with ϳ 0.2 in this case. The free energy barriers for the creation and annihilation of a droplet scales in 3D as ϳL , with ഛ ഛ 2.
Although there can be no dAT plane in the strict sense in the droplet-model, for an external magnetic field the system can be kept from equilibration on experimental time scales for parameters below a line that scales just like the dAT line. This phenomenon might mimic the effects of the replica symmetry breaking in the MP picture ͑see ͓107͔ for further details͒.

E. Replica-symmetric solution for fixed magnetization
This section serves to show that the mean-field version of the effective Hamiltonian ͑26͒ with random magnetic field and magnetization constraint exhibits replica symmetry breaking just as for the pure Edwards-Anderson model, and would therefore be a candidate to examine the validity of the MP or droplet picture in a realistic short-ranged spin glass model. Following Sherrington and Kirkpatrick ͑SK͒ ͓6͔, the mean-field model is given by where the round brackets ͑·,·͒ are generally used to denote sums over all pairs of different indices. This model differs from Eq. ͑26͒ by the long-range spin exchange. As the mean of i , H is generally nonzero this model will not exhibit a phase transition, which, however, is not a concern, as the number of ͑quasi-͒ground states will be the quantity of interest. Following the analysis of SK, we aim at finding the free energy, ground state overlap and magnetization constraint. Then we will use the de Almeida and Thouless approach ͓6͔ to show that the obtained solution is unstable in a certain parameter region, that lies below the so-called dAT-plane of stability. The type of instability that emerges is then wellknown to require the replica symmetry breaking solution of Parisi ͓6͔.
As the disorder is quenched ͑static on experimental time scales͒, one cannot average directly over disorder in the partition function as would be done for annealed disorder, but one must rather average the free energy density, f =−␤ ln Z using the "replica trick:" We form n identical copies of the system ͑the replicas͒ and the average is calculated for an integer n and a finite number of spins N. Then, using the general formula ln x = lim n→0 ͑x n −1͒ / n, ln Z is obtained from the analytic continuation of Z n for n → 0. Finally, we take the thermodynamic limit N → ϱ. Explicitly, Z n is given by where H SK ͓s i ␣ , n͔ is the sum of n independent and indentical spin Hamiltonians ͑29͒, averaged over the Gaussian disorder, with Greek indices now numbering the n replicas. Executing the average over the Gaussian distributions for K ij and i leads to coupling between spin-spin interactions of different replicas. As the mean-field approach means that the double sum over the site indices in Eq. ͑29͒ can be simplified into a square using ͑s i ␣ ͒ 2 = 1, one finds: where the prefactor e −n 2 ͑␤K͒ 2 /2+n͑␤h͒ 2 n/2 becomes irrelevant in the limit n → 0 and is subsequently dropped. As in the standard procedure, the square of the operator sum ͚ i s i ␣ s i ␤ is decoupled by introducing auxiliary operators q ␣␤ via a Hubbard-Stratonovitch ͑HS͒ transformation ͓6͔.
where the functional L͑q ␣␤ ͒ is and the configuration sum of exp͓L͑q ␣␤ ͔͒ now only goes over the n spins s ␣ in L͑q ␣␤ ͒, the HS transformation having made it possible to decouple the configuration sum over Nn spins in Eq. ͑31͒ into a N-fold product of n-spin sums. Assuming that the thermodynamic limit ͑N → ϱ͒ can be taken before n → 0, i.e., that the usual limiting process can be inverted, then Eq. ͑32͒ can be evaluated by the method of steepest descent, as the exponent is proportional to N. According to this method, the free energy per spin in the thermodynamic limit is the maximum of the q ␣␤ -dependent function in the exponent: ͑with f ª lim n→0 f n ͒ with the self-consistency condition: and the magnetization: where the average ͗͑·͒͘ L is defined as The mean-field approach has allowed a decoupling of the spins and a reduction of the problem to a single-site model with "Hamiltonian" L͓q ␣␤ ͔. For this new problem, the overlap parameter emerges naturally, albeit in a self-consistent manner. To push the calculation further some assumption for q ␣␤ has to be made. Naively, from the requirement that the result should be independent of the replica-indices, the most natural choice for q is to consider all identical overlaps between the replicas, q ␣␤ = q, which is the SK ansatz. Thus the double sum over the replicas ͚ ͑␣,␤͒ s ␣ s ␤ in Eq. ͑33͒ can be written as a square, keeping ͑s ␣ ͒ 2 = 1 in mind. Another Hubbard-Stratonovitch transformation with auxiliary vari-able z then decouples the square and yields an expression for the free energy density, which has to be evaluated selfconsistently and in which the limit n → 0 can be easily calculated: with A͑z͒ ª ␤ ͱ K 2 q + h 2 − ␤H. The overlap ͑35͒ and the magnetization constraint ͑36͒ can also be evaluated in the same way: A well-known problem with the SK ansatz for q ␣␤ is that it yields negative entropy for low temperatures and thus becomes unphysical. This is due to a fundamental technical problem with the replica trick: For the method of steepest descent to be valid, the SK solution must be a maximum of the exponent in Eq. ͑32͒ and must stay a maximum as the replica limit n → 0 is taken. But there is no unique way of choosing the zero-dimensional limit of the matrix q ␣␤ . The SK solution just corresponds to one possible choice for this limit. Thus the question arises whether the SK solution for the free energy is still a good, i.e., maximal choice in the replica limit.
To answer this question, one proceeds analogous to the de Almeida and Thouless procedure ͓137͔ to analyze the fluctuations around the SK solution, while taking the magnetization constraint into account ͑see Appendix B for details͒. Developing Eqs. ͑34͒ and ͑36͒ to second and first order, respectively, around q ␣␤ = q one finds that in the replica limit n → 0 there is an eigenvalue 2 of the matrix ‫ץ‬ 2 f / ‫ץ‬q ␣␤ ‫ץ‬q ␥␦ that can have both negative values and respects the constraint, yielding the condition which is violated for low enough T / K, H / K, and h / K. The plane in parameter space below which this happens is the so-called dAT surface. This instability is rectified by a much more involved ansatz for q ␣␤ that breaks the symmetry of the replicas and leads to the phenomena described in Sec. V A.

VI. EXPERIMENTAL INSIGHTS
Experimental creation and detection of the phenomena discussed in this paper poses a major problem and deserves a lot of creative thinking and separate publications. In this section we will just sketch what in our opinion are the most obvious ways of addressing these problems. In this section we do not address the questions concerning experimental realization of ultracold FB mixtures and composite fermions-these questions are discussed in Ref. ͓53͔.
The first question thus to be addressed is what are the best ways to create quenched disorder in a controlled way. Roth and Burnett ͓22͔ and the authors ͓21͔ have suggested that the use of pseudorandom disorder induced by noncommensurate optical lattices should work as well as the use of the genuine random lattices. Indeed the latter can be only ͑so far͒ achieved using speckle radiation, i.e., disorder correlation length of order of a few microns. If we work with systems of size of mm's, such disorder would definitely be enough to induce localization in 1D or 2D. Unfortunately, the size of the systems in question is typically of the order of hundred of microns, and that is one of the reasons why it is difficult to observe Anderson localization with BECs ͓27,28͔.
The analysis performed by us in this context implies clearly that it will be much easier to achieve the desired properties of the disorder using pseudodisordered, i.e., overlapped incommensurable optical lattices ͓29͔. One should also stress the equally promising look to the proposals formulated recently to use the optical tweezers techniques ͓38͔ and a random distribution of impurity atoms pinned in different lattice sites ͓33͔.
Another variety of problems is related to the detection of the quantum phases discussed in this paper. Below we list basic methods that have been already successfully applied to ultracold atomic gases in optical lattices.
͑1͒ Imaging of the atomic cloud after ballistic expansion. This ͑perhaps the most standard method͒ has been used in Ref. ͓20͔ to distinguish the bosonic SF phase from the MI phase. It allows for measurement of the quasimomentum distribution of atoms obtained after initial expansion caused by interactions ͓138͔, i.e., it detects first order coherence ͑interference pattern͒, present, for instance, in the SF phase.
͑2͒ Monitoring the density profile. Using phase contrast imaging ͓139͔ it is possible to perform a direct and nondestructive observation of the spatial distribution of the condensate in situ. This kind of measurement allows the direct observation of superfluidity ͓140͔ and can be therefore applied to characterize the fluid and superfluid phases.
͑3͒ Tilting or acceleration of the lattice. This method was used in Ref. ͓141͔ to detect the gap in the MI phase. It allows, in principle, to distinguish gapless from gapped phases, provided the continuum of low energy states can be achieved via tilting. Fluidity, superfluidity, and in general extended, nonlocalized excitations should allow one to detect Bloch oscillations ͓142͔.
͑4͒ Absorption of energy via modulation of the lattice. This method was also designed to detect the gap in the MI phase ͓124͔. Similarly as tilting, it provides a way of probing excitations in the systems.
͑5͒ Bragg spectroscopy. This method, one of the first proposed ͓143͔, is also a way of probing a certain kind of excitations in the system ͓144͔.
͑6͒ Cooper pair spectroscopy. This method is particularly useful to detect Fermi superfluids ͓145͔. Its theoretical aspects are discussed in Ref. ͓146͔ ͑9͒ Spatial quantum noise interferometry. The last, but not the least method discussed here allows for practically direct measurement of the density-density correlation and second order coherence. It has been proposed in Ref. ͓150͔ ͑see also ͓151͔͒, and used with great success to detect the bosonic Mott insulator ͓152͔ and the Fermi superfluid ͓153͔. It has capabilities of detecting and measuring relevant properties of various phases and structures ranging from supersolids, charge-density wave phases, and even Luttinger liquids in 1D. In the Fourier frequency and momentum domain it corresponds to measurements of the dynamic structure factor ͑for a discussion in the context of cold gases see Ref.

͓154͔͒.
There are of course more methods than the ones discussed above, but combined applications of those discussed should allow for clear detection and characterization of the quantum phases discussed in this paper.
Let us start this discussion with ideal Fermi gas, Fermi liquid, and metallic phases between Fermi glass and glassy Mott insulator. All of them are fluids, i.e., will respond consequently to perturbations. They are gapless and differ in this sense from the Fermi superfluids. All of them should lead to a nontrivial Fermi surface imaging in ballistic expansion. The difference between ideal and interacting phases is here rather quantitative, and as such can be measured. Quantities such as the effective mass can be recovered from the measurements. Influence of disorder on these phases will be seen as a gradual decrease of their "conducting" properties. A similar scenario is expected to take place with "dirty" superfluids; here measurements of the gap ͑using any of the excitation probing methods͒ should reveal a rapid gap decrease with the increasing disorder.
Disordered and glassy phases, such as Femi glass or glassy Mott insulator, are more difficult to detect. Obviously they will tend to give blurred images in the measurement of the first order coherence. Spatial noise interferometry should reveal some information about the glassy Mott insulator, especially in the region of parameters where it will incorporate domains of checkerboard phase. Although the glassy phases are gapless, the states forming the spectral quasicontinuum at low energies may be very difficult to achieve in simple excitation measurements, since they may lie far away one from another in the phase space. Probing of one of such states would thus allow one to study excitations accesible locally in the phase, which most presumably will be gapped. The character of excitations, and in particular their spectrum should, however, be a very sensitive function of how one excites them, and how one detects them ͑compare ͓32͔͒. On the other hand, the domain insulator phase should be visible by the "naked eye," and independent of detection of bosons and composites. Also, the noise interferometry should reveal information about the presence of the lattice, similarly as in the standard Mott insulator phase ͓152͔.
Finally, a separate problem concerns detection of the fermionic spin glass phase and its properties, as well as distinction between the possible adequacy of the Parisi versus droplet model. Repeated preparation of the system in the lattice with the same disorder ͑or even direct comparison and measurement of overlap of replicas ͓155͔͒ will shed some light on the latter problem. In many other aspects, response of the spin glass to excitations will be similar to that of Fermi glasses and glassy Mott insulator.

VII. CONCLUSIONS
Summarizing, we have studied atomic Fermi-Bose mixtures in optical lattices in the strong interacting limit, and in the presence of an inhomogeneous, or random on-site potential. We have derived the effective Hamiltonian describing the low temperature physics of the system, and shown that an inhomogeneous potential may be efficiently used to control the nature and strength of ͑boson mediated͒ interactions in the system. Using a random potential, one is able to control the system in such a way that its physics corresponds to a wide variety of quantum disordered systems. It is worth mentioning that the physics discussed in this paper is very much analogous to the one of Bose-Bose mixtures in the limit of hard core bosons ͑when both species exhibit strong intraspecies repulsion͒.
We end this section with a general comment on quantum complex systems. In our opinion quantum degenerate gases offer an absolutely unique possibility to study various models of physics of disorder systems, such as Bose and Fermi glasses, quantum percolating systems, "dirty" superfluids, domain and Mott insulators, quantum spin glasses, systems exhibiting localization-delocalization phenomena, etc. The summary of predictions of this paper is schematically shown in the following list of quantum phases, obtained for the case of only one type of sites, and i f =0, or i f = i b : ͑1͒ composites I-Fermi liquid, Fermi glass, quantum bond percolation; ͑2͒ composites II-ideal Fermi gas, Fermi glass, Fermi liquid, Mott insulator, fermionic spin glass; and ͑3͒ composites IĪ -domain insulator, "dirty" superfluid, Fermi glass, metallic phase, Mott insulator, fermionic spin glass. Additionally, for the case of lattices with different types of sites, physics of quantum site percolation will become relevant. Complex systems such as quantum cellular automata or neural networks can also be realized in this way. In fact, we and other authors have several times already stressed the fascinating possibility of using the ultracold lattice gases as quantum simulators of complex systems. But, the proposed systems go beyond just repeating what is known from the other kinds of physics; they allow one to create novel quantum phases and novel quantum behaviors.
͑2͒ H ef f does not couple states from different manifolds: P ␣ H ef f P ␤ = 0, ␣ ␤. ͑A2͒ ͑3͒ As the first two conditions still allow for an infinite number of unitary transformations ͑all UT are still possible, U being any unitary transformation acting only within the manifolds͒, the following additional condition is imposed: Expanding the first condition using the Baker-Hausdorff formula, one obtains: H ef f = H + ͓iS,H͔ + 1 2! ͓iS,͓iS,H͔͔ + 1 3! ͓iS,͓iS,͓iS,H͔͔͔ +¯. ͑A4͒ Making a power-series ansatz in S, S = JS 1 + J 2 S 2 + J 3 S 3 +¯, ͑A5͒ and employing H = H 0 + JH int one obtains from Eq. ͑A4͒ to second order

͑A6͒
This is a power series for H ef f , with its moments denoted by H ef f 1 , H ef f 2 ,..., where H ef f n generally depends on all S j with 1 ഛ j ഛ n. This allows for a systematic evaluation of the matrix elements ͗␣ , i͉iS j ͉␤ , j͘, and consequently delivers matrix expressions for the H ef f j and H ef f . To start with this, one considers the expansions ͑A6͒ and ͑A5͒ up to first order, i.e., H ef f = H 0 + JH ef f 1 and S = JS 1 . Using the second and third conditions, as well as P ␣ H 0 P ␤ = 0 and the expression for H ef f 1 in Eq. ͑A6͒, one finds: ͗␣,i͉iS 1 ͉␤, j͑͘E ␤j − E ␣i ͒ + ͗␣,i͉H int ͉␤, j͘ = 0 ͑A7͒ ⇒͗␣,i͉iS 1 ͉␤, j͘ = Ά ͗␣,i͉H int ͉␤, j͘ Thus the effective Hamiltonian within the ␣ manifold depends only on the interaction term and not on S 1 , i.e., ͗␣ , i͉H ef f 1 ͉␣ , j͘ = ͗␣ , i͉H int ͉␣ , j͘. A general result for any n is that ͗␣ , i͉H ef f n ͉␣ , j͘ is independent of S n , based on the third condition, and on the observation that S n enters the expression for H ef f n only in the commutator with H 0 , which is diagonal in the manifold index.
Thus when continuing to second order, the term ͓iS 2 , H 0 ͔ in the expression for H ef f 2 can be dropped. Of the two remaining terms defining H ef f 2 in Eq. ͑A6͒, the second one can be simplified by observing that according to Eq. ͑A7͒ the operator ͓iS 1 , H 0 ͔ is purely nondiagonal in the manifold index, with values opposite to those of the nondiagonal part of the interaction Hamiltonian. Thus 1 2 ͓iS 1 , ͓iS 1 , H 0 ͔͔ =− 1 2 ͓iS 1 , H int nd ͔. Now inserting the identity between the operators in the still untreated second term in H ef f 2 , ͓iS 1 , H int ͔, one sees that due to S 1 being nondiagonal in ␣, again only the nondiagonal part of H int can contribute: ͓iS 1 , H int ͔ = ͓iS 1 , H int nd ͔. Therefore one has

͑A10͒
where the identity operator has been inserted in the final expression for H ef f 2 in formula ͑A9͒, and then evaluated using formula ͑A8͒, which naturally leads one to define the operator Q ␣i as above. Note that this construction can be generalized to arbitrary orders in J in a straightforward manner, as detailed in ͓114͔.

APPENDIX B: STABILITY OF THE SK SOLUTION
The Taylor expansion of Eq. ͑34͒ around q ␣␤ = q yields the correction − 1 2 ␦f to the free energy, with
If the SK solution really corresponds to a maximum, this symmetric quadratic form must be positive definite. To check this, one calculates the eigenvalues of the 1 2 n͑n −1͒-dimensional matrix G ͑␣␤͒͑␥␦͒ in terms of its three distinct matrix elements: As becomes apparent from this, there are just three distinct classes of transformations that leave G ͑␣␤͒͑␥␦͒ invariant: those that permute no index acting on the P's; those that permute a single index acting on the Q's; and those that permute two indices acting on the R's. Thus there are only three eigenspaces to distinct eigenvalues and the linearly independent eigenvectors within each of these eigenspaces can naturally be chosen by considering a group of vectors that is invariant under the corresponding permutation transformation.

No index permutation
The ansatz for the eigenvalue is trivially given by ␦q ␣␤ = c, for all ͑␣␤͒, ͑B2͒ which is nondegenerate. With this ansatz, the eigenvalue equation is from which 1 can immediately be read off.

Permutation of a single index
Again, the ansatz for the eigenvectors is naturally given by ␦q ␣␤ = c for ␣ or ␤ = , ␦q ␣␤ = d for ␣,␤ .

͑B4͒
The ansatz still contains the previous case, as the requirement of the eigenvectors ͑B2͒ and ͑B4͒ being orthogonal still has to be fulfilled. This yields c = ͑1−n /2͒d and a degeneracy of n − 1 for the eigenvalue. The eigenvalue equation then becomes ͓P + ͑n − 4͒Q + ͑n − 3͒R − 1 Ј͔d = 0 ͑B5͒ from which 1 Ј is again immediately obvious.

͑B8͒
As de Almeida and Thouless report ͓137͔, a region in parameter-space where this limiting value took a negative value could not be found. Thus, calculating 2 explicitly, the relevant stability condition reads 2 = 1 ͑␤K͒ 2 − 1 ͱ 2 ͵ −ϱ ϱ dze −z 2 /2 sech 4 ͓A͑z͔͒ Ͼ 0. ͑B9͒ Solving the coupled equations ͑39͒-͑41͒ for 2 = 0 yields the dAT plane ͑cf. Fig. 11͒. Above it, the SK solution is still valid, and below it 2 takes negative value and the SK solution breaks down.

͑B11͒
Inserting any eigenvector of 2 into Eq. ͑B10͒ will, however, yield the desired result. As is clearly seen from Eq. ͑B11͒ ‫͗ץ‬s ͘ L / ‫ץ‬q ␣␤ evaluated at q ␣␤ = q is an eigenvector to 1 Ј, and thus ␦m = 0 is fulfilled.
The magnetization constraint is thereby compatible with the instability 2 Ͻ 0, and replica-symmetry breaking is expected to occur. From the calculations it has become clear that the presence of the random magnetic field would not change the occurrence of replica symmetry breaking, and all properties like infinitely degenerate ground-states and ultrametricity would be expected to occur in this model as well, provided the Mézard-Parisi approach can be applied to the finite range spin glass at all. FIG. 11. Almeida-Thouless plane in the reduced variables kT / K, h / K, and H / K. The constraint on the magnetization ͑41͒ has been effectively neglected for this, as it is a parameter in the experiment as well.