Difference between revisions of "ALPS 2 Tutorials:DMFT-02 Hybridization"

From ALPS
Jump to: navigation, search
(Tutorial 02: Hybridization Expansion CT-HYB)
(Tutorial 02: Hybridization Expansion CT-HYB)
Line 5: Line 5:
  
 
You can find the parameter files (called <tt>*tsqrt2</tt>) in the directory <tt>tutorials/dmft-02-hybridization</tt> in the examples.
 
You can find the parameter files (called <tt>*tsqrt2</tt>) in the directory <tt>tutorials/dmft-02-hybridization</tt> in the examples.
 +
Alternatively you can run the python script:
  
 +
import pyalps
 +
import numpy as np
 +
import matplotlib.pyplot as plt
 +
import pyalps.pyplot
 +
 +
#prepare the input parameters
 +
parms=[]
 +
for b in [6.,8.,10.,12.,14.,16.]:
 +
    parms.append(
 +
            {
 +
              'ANTIFERROMAGNET'    : 1,
 +
              'CONVERGED'          : 0.03,
 +
              'F'                  : 10,
 +
              'FLAVORS'            : 2,
 +
              'H'                  : 0,
 +
              'H_INIT'              : 0.2,
 +
              'MAX_IT'              : 10,
 +
              'MAX_TIME'            : 60,
 +
              'MU'                  : 0,
 +
              'N'                  : 1000,
 +
              'NMATSUBARA'          : 1000,
 +
              'N_FLIP'              : 0,
 +
              'N_MEAS'              : 10000,
 +
              'N_MOVE'              : 0,
 +
              'N_ORDER'            : 50,
 +
              'N_SHIFT'            : 0,
 +
              'OMEGA_LOOP'          : 1,
 +
              'OVERLAP'            : 0,
 +
              'SEED'                : 0,
 +
              'SITES'              : 1,
 +
              'SOLVER'              : 'Hybridization',
 +
              'SYMMETRIZATION'      : 0,
 +
              'TOLERANCE'          : 0.01,
 +
              'U'                  : 3,
 +
              't'                  : 0.707106781186547,
 +
              'SWEEPS'              : 100000000,
 +
              'THERMALIZATION'      : 1000,
 +
              'BETA'                : b,
 +
              'CHECKPOINT'          : 'solverdump_beta_'+str(b)+'.task1.out.h5',
 +
              'G0TAU_INPUT'        : 'G0_tau_input_beta_'+str(b)
 +
            }
 +
        )
 +
 +
#write the input file and run the simulation
 +
for p in parms:
 +
    input_file = pyalps.writeParameterFile('parm_beta_'+str(p['BETA']),p)
 +
    res = pyalps.runDMFT(input_file)
 +
 
After running these simulations compare the output to the Hirsch Fye results. To rerun a simulation, you can specify a starting solution by defining <tt>G0OMEGA_INPUT</tt>, e.g. copy <tt>G0omga_output</tt> to <tt>G0_omega_input</tt>, specify <tt>G0OMEGA_INPUT = G0_omega_input</tt> in the parameter file and rerun the code.  
 
After running these simulations compare the output to the Hirsch Fye results. To rerun a simulation, you can specify a starting solution by defining <tt>G0OMEGA_INPUT</tt>, e.g. copy <tt>G0omga_output</tt> to <tt>G0_omega_input</tt>, specify <tt>G0OMEGA_INPUT = G0_omega_input</tt> in the parameter file and rerun the code.  
 +
As in the Tutorial DMFT-01 Hirsch-Fye you can observe the transition to the antiferromagnetic phase by plotting the Green's functions using the script:
 +
 +
flavors=parms[0]['FLAVORS']
 +
listobs=[] 
 +
for f in range(0,flavors):
 +
    listobs.append('Green_'+str(f))
 +
 +
ll=pyalps.load.Hdf5Loader()
 +
data = ll.ReadMeasurementFromFile(pyalps.getResultFiles(pattern='parm_beta_*h5'), respath='/simulation/results/G_tau', measurements=listobs, verbose=True)
 +
for d in data:
 +
    for f in range(0,flavors):
 +
        d[f].x = d[f].x*d[f].props["BETA"]/float(d[f].props["N"])
 +
        plt.figure()
 +
        pyalps.pyplot.plot(d[f])
 +
        plt.xlabel(r'$\tau$')
 +
        plt.ylabel(r'$G(\tau)$')
 +
        plt.title('Hubbard model on the Bethe lattice')
 +
        plt.show()
  
 
You will notice that the results are relatively noisy. The reason for that is that the expansion order at such high temperatures is very small, which renders the measurement procedure inefficient. You can improve statistics by increasing the total run time (<tt>MAX_TIME</tt>) or by running it on more than one CPU. For running it with MPI, try <tt>mpirun -np procs dmft parameter_file</tt> or consult the man page of your mpi installation.
 
You will notice that the results are relatively noisy. The reason for that is that the expansion order at such high temperatures is very small, which renders the measurement procedure inefficient. You can improve statistics by increasing the total run time (<tt>MAX_TIME</tt>) or by running it on more than one CPU. For running it with MPI, try <tt>mpirun -np procs dmft parameter_file</tt> or consult the man page of your mpi installation.
 +
 +
If you want to check the convergence of your DMFT self-consistency, you can plot the Green's functions of different iterations using:
 +
 +
ll=pyalps.load.Hdf5Loader()
 +
for p in parms:
 +
    data = ll.ReadDMFTIterations(pyalps.getResultFiles(pattern='parm_beta_'+str(p['BETA'])+'.h5'), measurements=listobs, verbose=True)
 +
    grouped = pyalps.groupSets(pyalps.flatten(data), ['iteration'])
 +
    nd=[]
 +
    for group in grouped:
 +
        r = pyalps.DataSet()
 +
        r.y = np.array(group[0].y)
 +
        r.x = np.array([e*group[0].props['BETA']/float(group[0].props['N']) for e in group[0].x])
 +
        r.props = group[0].props
 +
        r.props['label'] = r.props['iteration']
 +
        nd.append( r )
 +
    plt.figure()
 +
    plt.xlabel(r'$\tau$')
 +
    plt.ylabel(r'$G(\tau)$')
 +
    plt.title(r'\beta = %.4s' %nd[0].props['BETA'])
 +
    pyalps.pyplot.plot(nd)
 +
    plt.legend()
 +
    plt.show()
  
 
= Tutorial 03: Interaction Expansion CT-INT =
 
= Tutorial 03: Interaction Expansion CT-INT =

Revision as of 22:18, 15 June 2010

Tutorial 02: Hybridization Expansion CT-HYB

We will now reproduce the same result with a continuous-time Quantum Monte Carlo code: the Hybridization or CT-HYB code by Werner et al. and Werner and Millis. The parameters are the same, apart from the command for the solver:

SOLVER=Hybridization

You can find the parameter files (called *tsqrt2) in the directory tutorials/dmft-02-hybridization in the examples. Alternatively you can run the python script:

import pyalps
import numpy as np
import matplotlib.pyplot as plt
import pyalps.pyplot

#prepare the input parameters
parms=[]
for b in [6.,8.,10.,12.,14.,16.]:
    parms.append(
            {
              'ANTIFERROMAGNET'     : 1,
              'CONVERGED'           : 0.03,
              'F'                   : 10,
              'FLAVORS'             : 2,
              'H'                   : 0,
              'H_INIT'              : 0.2,
              'MAX_IT'              : 10,
              'MAX_TIME'            : 60,
              'MU'                  : 0,
              'N'                   : 1000,
              'NMATSUBARA'          : 1000,
              'N_FLIP'              : 0,
              'N_MEAS'              : 10000,
              'N_MOVE'              : 0,
              'N_ORDER'             : 50,
              'N_SHIFT'             : 0,
              'OMEGA_LOOP'          : 1,
              'OVERLAP'             : 0,
              'SEED'                : 0,
              'SITES'               : 1,
              'SOLVER'              : 'Hybridization',
              'SYMMETRIZATION'      : 0,
              'TOLERANCE'           : 0.01,
              'U'                   : 3,
              't'                   : 0.707106781186547,
              'SWEEPS'              : 100000000,
              'THERMALIZATION'      : 1000,
              'BETA'                : b,
              'CHECKPOINT'          : 'solverdump_beta_'+str(b)+'.task1.out.h5',
              'G0TAU_INPUT'         : 'G0_tau_input_beta_'+str(b)
            }
        )
#write the input file and run the simulation
for p in parms:
    input_file = pyalps.writeParameterFile('parm_beta_'+str(p['BETA']),p)
    res = pyalps.runDMFT(input_file)

After running these simulations compare the output to the Hirsch Fye results. To rerun a simulation, you can specify a starting solution by defining G0OMEGA_INPUT, e.g. copy G0omga_output to G0_omega_input, specify G0OMEGA_INPUT = G0_omega_input in the parameter file and rerun the code. As in the Tutorial DMFT-01 Hirsch-Fye you can observe the transition to the antiferromagnetic phase by plotting the Green's functions using the script:

flavors=parms[0]['FLAVORS']
listobs=[]   
for f in range(0,flavors):
    listobs.append('Green_'+str(f))

ll=pyalps.load.Hdf5Loader()
data = ll.ReadMeasurementFromFile(pyalps.getResultFiles(pattern='parm_beta_*h5'), respath='/simulation/results/G_tau', measurements=listobs, verbose=True)
for d in data:
    for f in range(0,flavors):
        d[f].x = d[f].x*d[f].props["BETA"]/float(d[f].props["N"])
        plt.figure()
        pyalps.pyplot.plot(d[f])
        plt.xlabel(r'$\tau$')
        plt.ylabel(r'$G(\tau)$')
        plt.title('Hubbard model on the Bethe lattice')
        plt.show()

You will notice that the results are relatively noisy. The reason for that is that the expansion order at such high temperatures is very small, which renders the measurement procedure inefficient. You can improve statistics by increasing the total run time (MAX_TIME) or by running it on more than one CPU. For running it with MPI, try mpirun -np procs dmft parameter_file or consult the man page of your mpi installation.

If you want to check the convergence of your DMFT self-consistency, you can plot the Green's functions of different iterations using:

ll=pyalps.load.Hdf5Loader()
for p in parms:
    data = ll.ReadDMFTIterations(pyalps.getResultFiles(pattern='parm_beta_'+str(p['BETA'])+'.h5'), measurements=listobs, verbose=True)
    grouped = pyalps.groupSets(pyalps.flatten(data), ['iteration'])
    nd=[]
    for group in grouped:
        r = pyalps.DataSet()
        r.y = np.array(group[0].y)
        r.x = np.array([e*group[0].props['BETA']/float(group[0].props['N']) for e in group[0].x])
        r.props = group[0].props
        r.props['label'] = r.props['iteration']
        nd.append( r )
    plt.figure()
    plt.xlabel(r'$\tau$')
    plt.ylabel(r'$G(\tau)$')
    plt.title(r'\beta = %.4s' %nd[0].props['BETA'])
    pyalps.pyplot.plot(nd)
    plt.legend()
    plt.show()

Tutorial 03: Interaction Expansion CT-INT

It is also instructive to run these calculations with a CT-INT code. This code performs an expansion in the interaction (instead of the hybridization). The corresponding parameter files are very similar, you can find them in the directory tutorials/dmft-03-interaction.

Tutorial by Emanuel - Please don't hesitate to ask!