#######################
Defining lattice models 
#######################

This section explains how to define lattice models through python calls.

A simple example
----------------

The following simple example illustrates how to define a Hubbard model on the square lattice in two dimensions::

    import pyqcm

    CM = pyqcm.cluster_model(4, 0, 'clus', ((4, 3, 2, 1),))
    clus = pyqcm.cluster(CM, ((0, 0, 0), (1, 0, 0), (0, 1, 0), (1, 1, 0)))
    model = pyqcm.lattice_model('2x2_C2', clus, ((2, 0, 0), (0, 2, 0)))
    model.interaction_operator('U')
    model.hopping_operator('t', (1, 0, 0), -1)
    model.hopping_operator('t', (0, 1, 0), -1)
    model.hopping_operator('t2', (1, 1, 0), -1)
    model.hopping_operator('t2', (-1, 1, 0), -1)
    model.hopping_operator('t3', (2, 0, 0), -1)
    model.hopping_operator('t3', (0, 2, 0), -1)
    model.anomalous_operator('D', (1, 0, 0), 1)
    model.anomalous_operator('D', (0, 1, 0), -1)
    model.density_wave('M', 'Z', (1, 1, 0))

There is a call to ``cluster_model()`` to define a :math:`2\times2` plaquette with a :math:`C_2` symmetry (rotation by 180 degrees). Then another object, of type ``cluster``, is defined based on this model. This object has geometrical meaning, whereas the ``cluster_model`` does not. The ``cluster`` constructor takes the following arguments:

#. The ``cluster_model`` object it is based on. It is possible to associate more than one cluster model to a cluster, if one wants to use the *subbath* concept is dynamical mean field theory.
#. An array of integer positions of the different cluster sites. This is where the geometry of the cluster appears.
#. The base position of the cluster within the super unit cell ( ``[0,0,0]`` by default). This is just added for convenience when complex super unit cells are built out of the same cluster, shifted in position.

Once the different clusters (or the single cluster) have been defined, one can construct an object of type ``lattice_model`` that contains the different clusters arranged in a super-lattice. The constructor takes the following arguments:

#. The name of the model, for reporting purposes
#. The cluster (or sequence of clusters) forming the repeated unit.
#. The :math:`d` vectors defining the superlattice vectors. :math:`d` is the dimension of the model.
#. The :math:`d` vectors defining the lattice vectors of the model. In the above example this is omitted because the default value is used: ``((1,0,0), (0,1,0))``.

Note that all vectors used in these specifications have three components (even for lower dimensional models) and have **integer** components. They are in fact coefficients of basis vectors generating a working Bravais lattice and are therefore integers by definition (not all vectors of this working Bravais lattice need belong to the actual Bravais lattice of the model; the case of the graphene lattice is illustrated below). They can be tuples, lists or numpy arrays.
We recommend using tuples for simplicity.

Once the geometry of the model is defined, operators can be added to the model via various functions. The Hubbard on-site interaction is added with a call to the ``lattice_model`` member function  ``interaction_operator('U')``, whose sole argument in this case is the name we choose for that operator (more arguments are needed for other types of interactions, multi-band models, etc; see the detailed documents in the reference section on functions).

Nearest-neighbor hopping is defined with a call to the member function ``hopping_operator()``, which has three mandatory arguments:

#. The name of the operator
#. The link on which the operator is defined
#. The amplitude of the operator on that link.

Optional keyword arguments are needed for multi-band models, spin-flip operators, etc. In the above example this function is called twice with the same name but different links (along (1,0,0) and (0,1,0) respectively). The matrix elements generated by the two calls are simply added to the list associated with this operator. Similar calls are performed for the second-neighbor hopping ``t2`` and the third-neighbor hopping ``t3``.

A d-wave pairing operator is defined via a call to ``anomalous_operator()`` , which takes the same arguments as ``hopping_operator()``, i.e., name, link and amplitude. Note that two calls are made with the same name ``D``, one in the :math:`x` direction, the other one in the :math:`y` direction, with amplitudes 1 and -1, in accordance with the d-wave character of the operator.


Finally, a density wave corresponding to :math:`(\pi,\pi)` antiferromagnetism is added, with a call to ``density_wave()``, which has 3 mandatory arguments:

#. The name of the operator (here 'M')
#. The type of density-wave (here 'Z' for a spin density wave in the :math:`z` component of the spin). Other possibilities are:
    * 'N' : charge density wave 
    * 'X' : spin density wave in the :math:`x` component of the spin
    * 'spin' : same as `Z`
    * 'singlet' : singlet pairing
    * 'dx' : triplet pairing in the :math:`x` direction using the :math:`\mathbf{d}` vector formalism.
    * 'dy' : same, in the :math:`y` direction
    * 'dz' : same, in the :math:`z` direction (this one does not require Nambu doubling)
#. The wavevector :math:`\mathbf{Q}` of the density wave. It has real components (in multiples of :math:`\pi`) and it must be commensurate with the super unit cell within some small numerical tolerance.

Additional keyword arguments to ``density_wave()`` include the link on which the density wave is defined (for bond-density waves), lattice orbitals involved (for multi-band models), additional phases and amplitudes, etc. Again, see the reference section for details.


A more complex example
----------------------

The following example defines a model on the graphene lattice using two clusters within the super unit cell and the graphene lattice, as illustrated below:

.. figure:: h8.png
    :align: center
    :height: 300px

    Figure 3

A possible set of function calls to define the Hubbard model on this system is::

    from cluster_h4_6b_C3 import CM
    import pyqcm

    C1 = pyqcm.cluster(CM, [[0, 0, 0], [-1, 0, 0], [1, 1, 0], [0, -1, 0]], [-1, 0, 0])
    C2 = pyqcm.cluster(CM, [[0, 0, 0], [1, 0, 0], [-1, -1, 0], [0, 1, 0]], [1, 0, 0])
    model = pyqcm.lattice_model('h4_6b_C3', (C1, C2), [[2, -2, 0],[2, 4, 0]], [[1, -1, 0], [2, 1, 0]])
    model.set_basis([[1, 0, 0], [-0.5, 0.866025403784438, 0], [0, 0, 1]])

    model.interaction_operator('U')
    model.hopping_operator('t', [1,0,0], -1, orbitals=(2,1))
    model.hopping_operator('t', [0,1,0], -1, orbitals=(2,1))
    model.hopping_operator('t', [-1,-1,0], -1, orbitals=(2,1))

The first statement, ``from cluster_h4_6b_C3 import CM``, imports a cluster definition file, for instance the one associated with Fig. 2 above, in which the ``cluster_model`` object named ``CM`` was defined. 
Two copies of the clusters are added to the super unit cell. The positions associated with the two copies are different, but the cluster model is the same, which means that only one copy of the Hilbert space operators and bases necessary for the exact diagonalization will be constructed. The origin has been placed exactly between the two clusters. The positions in each cluster are defined relative to the base position of each cluster.

The integer positions are defined in terms of the basis defined by the call to ``set_basis()``. The argument of that function is a set of real-valued vectors defining the basis vectors of the working Bravais lattice.
On Fig. 3, the first two of these vectors are :math:`\mathbf{e}_1` and :math:`\mathbf{e}_2`.

The call to ``lattice_model()`` defines both the superlattice vectors :math:`\mathbf{E}_1` and :math:`\mathbf{E}_2` (second argument) and the basis vectors of the model's physical Bravais lattice (third argument).
The lattice basis vectors only serve to attribute orbital labels to the different sites. In `qcm`, each degree of freedom of a given spin (i.e. each orbital) must have its own site on the working Bravais lattice. The basis vectors of the physical Bravais lattice then define orbital labels (from 1 to :math:`N_\mathrm{band}`) attributed to each site. The order in which orbitals are labelled depends on the order in which the sites appear. Given the above definitions and Fig. 3, the A sublattice of graphene corresponds to orbital 2 and the B sublattice to orbital 1.

Given that the model has two bands, the definition of the hopping operator ``t`` must contain orbital information: the keywords ``orbitals`` (a tuple) is used to specify the orbital numbers associated with the two sites separated by the bond vector (link) given in argument. Internally, a loop is done over all sites of the super unit cell; the bond vector is used to identify a second site; if that site exists and if the orbitals associated with the two sites agree with ``orbitals``, then a matrix element is added to the operator. In the above example, three calls are needed because of the three directions (bonds). 

The greatest risk in such calls is to mislabel the orbitals. In order to check that operators were defined properly, a 
submodule called ``draw_operator`` is provided, that can draw on the screen a schematic view of each operator defined on the lattice or in a particular cluster model.

Other examples
--------------

The distribution contains a folder (`notebooks`) that contains many examples of models and codes. New users are encouraged to study a few of these models and to consult the reference section for more detailed information about model building.

Class lattice_model
-------------------

.. autoclass:: pyqcm.lattice_model 
    :members: