From a882ffa5759e5a781600b223960df3d90ab57e56 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Malte=20Sch=C3=BCler?= Date: Mon, 8 Jul 2019 09:55:22 +0200 Subject: [PATCH] work in NiO tutorial --- doc/CMakeLists.txt | 2 +- doc/guide/conv_vasp.rst | 9 +- doc/guide/plovasp.rst | 21 ++- doc/tutorials.rst | 9 ++ doc/tutorials/images_scripts/INCAR | 18 +++ doc/tutorials/images_scripts/INCAR.rst | 7 + doc/tutorials/images_scripts/KPOINTS | 4 + doc/tutorials/images_scripts/KPOINTS.rst | 7 + .../images_scripts/NiO_local_lattice_GF.py | 90 +++++++++++ doc/tutorials/images_scripts/POSCAR | 10 ++ doc/tutorials/images_scripts/POSCAR.rst | 7 + doc/tutorials/images_scripts/converter.py | 3 + doc/tutorials/images_scripts/converter.py.rst | 10 ++ doc/tutorials/images_scripts/nio.py | 149 ++++++++++++++++++ doc/tutorials/images_scripts/nio.py.rst | 9 ++ doc/tutorials/images_scripts/plo.cfg | 27 ++++ doc/tutorials/images_scripts/plo.cfg.rst | 9 ++ doc/tutorials/nio_csc.rst | 59 +++++++ python/converters/vasp_converter.py | 2 +- python/sumk_dft.py | 8 +- test/plovasp/converter/lunio3.out.h5 | Bin 336304 -> 372392 bytes test/plovasp/converter/pg_output.out.h5 | Bin 221512 -> 282480 bytes 22 files changed, 445 insertions(+), 15 deletions(-) create mode 100644 doc/tutorials/images_scripts/INCAR create mode 100644 doc/tutorials/images_scripts/INCAR.rst create mode 100644 doc/tutorials/images_scripts/KPOINTS create mode 100644 doc/tutorials/images_scripts/KPOINTS.rst create mode 100644 doc/tutorials/images_scripts/NiO_local_lattice_GF.py create mode 100644 doc/tutorials/images_scripts/POSCAR create mode 100644 doc/tutorials/images_scripts/POSCAR.rst create mode 100644 doc/tutorials/images_scripts/converter.py create mode 100644 doc/tutorials/images_scripts/converter.py.rst create mode 100644 doc/tutorials/images_scripts/nio.py create mode 100644 doc/tutorials/images_scripts/nio.py.rst create mode 100644 doc/tutorials/images_scripts/plo.cfg create mode 100644 doc/tutorials/images_scripts/plo.cfg.rst create mode 100644 doc/tutorials/nio_csc.rst diff --git a/doc/CMakeLists.txt b/doc/CMakeLists.txt index af1d88d8..0e137f0e 100644 --- a/doc/CMakeLists.txt +++ b/doc/CMakeLists.txt @@ -18,6 +18,6 @@ add_custom_target(doc_sphinx ALL DEPENDS ${sphinx_top} ${CMAKE_CURRENT_BINARY_DI # --------------------------------- install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/html/ COMPONENT documentation DESTINATION share/doc/triqs_dft_tools FILES_MATCHING - REGEX "\\.(html|pdf|png|gif|jpg|js|xsl|css|py|txt|inv|bib)$" + REGEX "\\.(html|pdf|png|gif|jpg|js|xsl|css|py|txt|inv|bib|cfg)$" PATTERN "_*" ) diff --git a/doc/guide/conv_vasp.rst b/doc/guide/conv_vasp.rst index fb80eac7..9064cfd9 100644 --- a/doc/guide/conv_vasp.rst +++ b/doc/guide/conv_vasp.rst @@ -5,9 +5,7 @@ Interface with VASP =================== .. warning:: - The VASP interface is in the alpha-version and the VASP part of it is not - yet publicly released. The documentation may, thus, be subject to changes - before the final release. + The VASP interface is in the alpha-version. The documentation may, thus, be subject to changes before the final release. *Limitations of the alpha-version:* @@ -26,7 +24,8 @@ in the :ref:`PLOVasp User's Guide `. Here, a quick-start guide is prese The VASP interface relies on new options introduced since version 5.4.x. In particular, a new INCAR-option `LOCPROJ` -and new `LORBIT` modes 13 and 14 have been added. +and new `LORBIT` modes 13 and 14 have been added as well as the new ICHARG +mode 5 for charge self-consistent calculations Option `LOCPROJ` selects a set of localized projectors that will be written to file `LOCPROJ` after a successful VASP run. @@ -41,7 +40,7 @@ with the indices corresponding to the site position in the POSCAR file; `` chooses a particular type of the local basis function. The recommended projector type is `Pr 2`. The formalism for this type of projectors is presented in -`M. Schüler et al. 2018 J. Phys.: Condens. Matter 30 475901 `_. +`M. Schüler et al. 2018 J. Phys.: Condens. Matter 30 475901 `_. For details on `LOCPROJ` also have a look in the `VASP wiki `_ The allowed labels of the local states defined in terms of cubic harmonics are: diff --git a/doc/guide/plovasp.rst b/doc/guide/plovasp.rst index 99f5e1d9..b1e64fe8 100644 --- a/doc/guide/plovasp.rst +++ b/doc/guide/plovasp.rst @@ -55,21 +55,25 @@ A PLOVasp input file can contain three types of sections: Section [General] """"""""""""""""" -The entire section is optional and it contains three parameters: +The entire section is optional and it contains four parameters: * **BASENAME** (string): provides a base name for output files. Default filenames are :file:`vasp.*`. * **DOSMESH** ([float float] integer): if this parameter is given, the projected density of states for each projected orbital will be - evaluated and stored to files :file:`pdos_.dat`, where `n` is the - orbital index. The energy - mesh is defined by three numbers: `EMIN` `EMAX` `NPOINTS`. The first two + evaluated and stored to files :file:`pdos__.dat`, where `s` is the + shell index and `n` the ion index. The energy mesh is defined by three + numbers: `EMIN` `EMAX` `NPOINTS`. The first two can be omitted in which case they are taken to be equal to the projector energy window. **Important note**: at the moment this option works only if the tetrahedron integration method (`ISMEAR = -4` or `-5`) is used in VASP to produce `LOCPROJ`. * **EFERMI** (float): provides the Fermi level. This value overrides the one extracted from VASP output files. +* **HK** (True/False): If True, the projectors are applied the the Kohn-Sham + eigenvalues which results in a Hamitlonian H(k) in orbital basis. The H(k) + is written for each group to a file :file:`Basename.hk`. It is recommended + to also set `COMPLEMENT = True` (see below). Default is False. There are no required parameters in this section. @@ -94,6 +98,8 @@ given by `LSHELL` must be present in the LOCPROJ file. There are additional optional parameters that allow one to transform the local states: +* **CORR** (True/False): Determines if shell is correlated or not. At least one + shell has to be correlated. Default is True. * **TRANSFORM** (matrix): local transformation matrix applied to all states in the projector shell. The matrix is defined by a (multiline) block of floats, with each line corresponding to a row. The number of columns @@ -131,6 +137,13 @@ Optional group parameters: condition will be enforced on each site separately but the Wannier functions on different sites will not be orthogonal. If `NORMION = False`, the Wannier functions on different sites included in the group will be orthogonal to each other. +* **BANDS** (int int): the energy window specified by two ints: band index of + lowest band and band index of highest band. Using this overrides the selection + in `EWINDOW`. +* **COMPLEMENT** (True/False). If True, the orthogonal complement is calculated + resulting in unitary (quadratic) projectors, i.e., the same number of orbitals + as bands. It is required to have an equal number of bands in the energy window + at each k-point. Default is False. .. _transformation_file: diff --git a/doc/tutorials.rst b/doc/tutorials.rst index d02cdeb6..dd8c7a15 100644 --- a/doc/tutorials.rst +++ b/doc/tutorials.rst @@ -23,3 +23,12 @@ Full charge self consistency with Wien2k: :math:`\gamma`-Ce tutorials/ce-gamma-fscs_wien2k +A full example with VASP: NiO +----------------------------------------------------------- + +.. toctree:: + :maxdepth: 2 + + + tutorials/nio_csc + diff --git a/doc/tutorials/images_scripts/INCAR b/doc/tutorials/images_scripts/INCAR new file mode 100644 index 00000000..45dc19f2 --- /dev/null +++ b/doc/tutorials/images_scripts/INCAR @@ -0,0 +1,18 @@ +System = NiO + +ISMEAR = -5 + +# the energy window to optimize projector channels +EMIN = -3 +EMAX = 7 + +NBANDS = 16 +LMAXMIX = 6 + +# switch off all symmetries +ISYM = -1 + +# project to Ni d and O p states +LORBIT = 14 +LOCPROJ = 1 : d : Pr 1 +LOCPROJ = 2 : p : Pr 1 diff --git a/doc/tutorials/images_scripts/INCAR.rst b/doc/tutorials/images_scripts/INCAR.rst new file mode 100644 index 00000000..0a5cc0ac --- /dev/null +++ b/doc/tutorials/images_scripts/INCAR.rst @@ -0,0 +1,7 @@ +.. _INCAR: + +INCAR +----------- + +.. literalinclude:: INCAR + :language: bash diff --git a/doc/tutorials/images_scripts/KPOINTS b/doc/tutorials/images_scripts/KPOINTS new file mode 100644 index 00000000..7099c814 --- /dev/null +++ b/doc/tutorials/images_scripts/KPOINTS @@ -0,0 +1,4 @@ +Automatically generated mesh + 0 +Gamma + 6 6 6 \ No newline at end of file diff --git a/doc/tutorials/images_scripts/KPOINTS.rst b/doc/tutorials/images_scripts/KPOINTS.rst new file mode 100644 index 00000000..93aa5f3b --- /dev/null +++ b/doc/tutorials/images_scripts/KPOINTS.rst @@ -0,0 +1,7 @@ +.. _KPOINTS: + +KPOINTS +----------- + +.. literalinclude:: KPOINTS + :language: bash diff --git a/doc/tutorials/images_scripts/NiO_local_lattice_GF.py b/doc/tutorials/images_scripts/NiO_local_lattice_GF.py new file mode 100644 index 00000000..02dbe051 --- /dev/null +++ b/doc/tutorials/images_scripts/NiO_local_lattice_GF.py @@ -0,0 +1,90 @@ +from itertools import * +import numpy as np +import pytriqs.utility.mpi as mpi +from pytriqs.archive import * +from pytriqs.gf import * +from triqs_dft_tools.sumk_dft import * +from triqs_dft_tools.sumk_dft_tools import * +from pytriqs.operators.util.hamiltonians import * +from pytriqs.operators.util.U_matrix import * +from triqs_cthyb import * +import warnings +warnings.filterwarnings("ignore", category=FutureWarning) + +filename = 'vasp' +SK = SumkDFT(hdf_file = filename+'.h5', use_dft_blocks = False) + +beta = 5.0 + +# We analyze the block structure of the Hamiltonian +Sigma = SK.block_structure.create_gf(beta=beta) + +SK.put_Sigma([Sigma]) +G = SK.extract_G_loc() +SK.analyse_block_structure_from_gf(G, threshold = 1e-3) + + +# Setup CTQMC Solver +n_orb = SK.corr_shells[0]['dim'] +spin_names = ['up','down'] +orb_names = [i for i in range(0,n_orb)] + + +# Print some information on the master node +iteration_offset = 0 +Sigma_iw = 0 +if mpi.is_master_node(): + ar = HDFArchive(filename+'.h5','a') + if not 'DMFT_results' in ar: ar.create_group('DMFT_results') + if not 'Iterations' in ar['DMFT_results']: ar['DMFT_results'].create_group('Iterations') + if 'iteration_count' in ar['DMFT_results']: + iteration_offset = ar['DMFT_results']['iteration_count']+1 + print('offset',iteration_offset) + Sigma_iw = ar['DMFT_results']['Iterations']['Sigma_it'+str(iteration_offset-1)] + SK.dc_imp = ar['DMFT_results']['Iterations']['dc_imp'+str(iteration_offset-1)] + SK.dc_energ = ar['DMFT_results']['Iterations']['dc_energ'+str(iteration_offset-1)] + SK.chemical_potential = ar['DMFT_results']['Iterations']['chemical_potential'+str(iteration_offset-1)].real + +iteration_offset = mpi.bcast(iteration_offset) +Sigma_iw = mpi.bcast(Sigma_iw) +SK.dc_imp = mpi.bcast(SK.dc_imp) +SK.dc_energ = mpi.bcast(SK.dc_energ) +SK.chemical_potential = mpi.bcast(SK.chemical_potential) + + +SK.put_Sigma(Sigma_imp = [Sigma_iw]) + +ikarray = numpy.array(range(SK.n_k)) + +# set up the orbitally resolved local lattice greens function: +n_orbs = SK.proj_mat_csc.shape[2] +spn = SK.spin_block_names[SK.SO] +mesh = Sigma_iw.mesh +block_structure = [range(n_orbs) for sp in spn] +gf_struct = [(spn[isp], block_structure[isp]) + for isp in range(SK.n_spin_blocks[SK.SO])] +block_ind_list = [block for block, inner in gf_struct] +glist = lambda: [GfImFreq(indices=inner, mesh=mesh) + for block, inner in gf_struct] + +G_latt_orb = BlockGf(name_list=block_ind_list, + block_list=glist(), make_copies=False) +G_latt_orb.zero() + +for ik in mpi.slice_array(ikarray): + G_latt_KS = SK.lattice_gf(ik=ik, beta=beta) + G_latt_KS *= SK.bz_weights[ik] + for bname, gf in G_latt_orb: + add_g_ik = gf.copy() + add_g_ik.zero() + add_g_ik << SK.downfold(ik, 0, bname, G_latt_KS[bname], gf, shells='csc', ir=None) + gf << gf + add_g_ik + +G_latt_orb << mpi.all_reduce( + mpi.world, G_latt_orb, lambda x, y: x + y) + +mpi.barrier() + +if mpi.is_master_node(): + ar['DMFT_results']['Iterations']['G_latt_orb_it'+str(iteration_offset-1)] = G_latt_orb +if mpi.is_master_node(): del ar diff --git a/doc/tutorials/images_scripts/POSCAR b/doc/tutorials/images_scripts/POSCAR new file mode 100644 index 00000000..65f73f29 --- /dev/null +++ b/doc/tutorials/images_scripts/POSCAR @@ -0,0 +1,10 @@ +NiO +4.17 #taken from 9x9x9 with sigma=0.2 ismear=2 + +0.0000000000 +0.5000000000 +0.5000000000 + +0.5000000000 +0.0000000000 +0.5000000000 + +0.5000000000 +0.5000000000 +0.0000000000 +Ni O + 1 1 +Direct + +0.0000000000 +0.0000000000 +0.0000000000 + +0.5000000000 +0.5000000000 +0.5000000000 diff --git a/doc/tutorials/images_scripts/POSCAR.rst b/doc/tutorials/images_scripts/POSCAR.rst new file mode 100644 index 00000000..35153ff5 --- /dev/null +++ b/doc/tutorials/images_scripts/POSCAR.rst @@ -0,0 +1,7 @@ +.. _POSCAR: + +POSCAR +----------- + +.. literalinclude:: POSCAR + :language: bash diff --git a/doc/tutorials/images_scripts/converter.py b/doc/tutorials/images_scripts/converter.py new file mode 100644 index 00000000..8bdabb9b --- /dev/null +++ b/doc/tutorials/images_scripts/converter.py @@ -0,0 +1,3 @@ +from triqs_dft_tools.converters.vasp_converter import * +Converter = VaspConverter(filename = 'nio') +Converter.convert_dft_input() diff --git a/doc/tutorials/images_scripts/converter.py.rst b/doc/tutorials/images_scripts/converter.py.rst new file mode 100644 index 00000000..4b1ca3a8 --- /dev/null +++ b/doc/tutorials/images_scripts/converter.py.rst @@ -0,0 +1,10 @@ +.. _converter.py: + +converter.py +------------------- + +Download :download:`converter.py <./converter.py>`. + + +.. literalinclude:: converter.py + :language: python diff --git a/doc/tutorials/images_scripts/nio.py b/doc/tutorials/images_scripts/nio.py new file mode 100644 index 00000000..1d238840 --- /dev/null +++ b/doc/tutorials/images_scripts/nio.py @@ -0,0 +1,149 @@ +from itertools import * +import numpy as np +import pytriqs.utility.mpi as mpi +from pytriqs.archive import * +from pytriqs.gf import * +from triqs_dft_tools.sumk_dft import * +from triqs_dft_tools.sumk_dft_tools import * +from pytriqs.operators.util.hamiltonians import * +from pytriqs.operators.util.U_matrix import * +from triqs_cthyb import * + +filename = 'vasp' + +SK = SumkDFT(hdf_file = filename+'.h5', use_dft_blocks = False) + +beta = 5.0 + +Sigma = SK.block_structure.create_gf(beta=beta) +SK.put_Sigma([Sigma]) +G = SK.extract_G_loc() +SK.analyse_block_structure_from_gf(G, threshold = 1e-3) +for i_sh in range(len(SK.deg_shells)): + num_block_deg_orbs = len(SK.deg_shells[i_sh]) + mpi.report('found {0:d} blocks of degenerate orbitals in shell {1:d}'.format(num_block_deg_orbs, i_sh)) + for iblock in range(num_block_deg_orbs): + mpi.report('block {0:d} consists of orbitals:'.format(iblock)) + for keys in SK.deg_shells[i_sh][iblock].keys(): + mpi.report(' '+keys) + +# Setup CTQMC Solver + +n_orb = SK.corr_shells[0]['dim'] +spin_names = ['up','down'] +orb_names = [i for i in range(0,n_orb)] +orb_hyb = False + +#gf_struct = set_operator_structure(spin_names, orb_names, orb_hyb) +gf_struct = SK.gf_struct_solver[0] +mpi.report('Sumk to Solver: %s'%SK.sumk_to_solver) +mpi.report('GF struct sumk: %s'%SK.gf_struct_sumk) +mpi.report('GF struct solver: %s'%SK.gf_struct_solver) + +S = Solver(beta=beta, gf_struct=gf_struct) + +# Construct the Hamiltonian and save it in Hamiltonian_store.txt +H = Operator() +U = 8.0 +J = 1.0 + + +U_sph = U_matrix(l=2, U_int=U, J_hund=J) +U_cubic = transform_U_matrix(U_sph, spherical_to_cubic(l=2, convention='')) +Umat, Upmat = reduce_4index_to_2index(U_cubic) + +H = h_int_density(spin_names, orb_names, map_operator_structure=SK.sumk_to_solver[0], U=Umat, Uprime=Upmat) + +# Print some information on the master node +mpi.report('Greens function structure is: %s '%gf_struct) +mpi.report('U Matrix set to:\n%s'%Umat) +mpi.report('Up Matrix set to:\n%s'%Upmat) + +# Parameters for the CTQMC Solver +p = {} +p["max_time"] = -1 +p["random_name"] = "" +p["random_seed"] = 123 * mpi.rank + 567 +p["length_cycle"] = 100 +p["n_warmup_cycles"] = 20000 +p["n_cycles"] = 200000 +p["fit_max_moment"] = 4 +p["fit_min_n"] = 30 +p["fit_max_n"] = 50 +p["perform_tail_fit"] = True + +# Double Counting: 0 FLL, 1 Held, 2 AMF +DC_type = 0 + +# Prepare hdf file and and check for previous iterations +n_iterations = 10 + +iteration_offset = 0 +if mpi.is_master_node(): + ar = HDFArchive(filename+'.h5','a') + if not 'DMFT_results' in ar: ar.create_group('DMFT_results') + if not 'Iterations' in ar['DMFT_results']: ar['DMFT_results'].create_group('Iterations') + if 'iteration_count' in ar['DMFT_results']: + iteration_offset = ar['DMFT_results']['iteration_count']+1 + S.Sigma_iw = ar['DMFT_results']['Iterations']['Sigma_it'+str(iteration_offset-1)] + SK.dc_imp = ar['DMFT_results']['Iterations']['dc_imp'+str(iteration_offset-1)] + SK.dc_energ = ar['DMFT_results']['Iterations']['dc_energ'+str(iteration_offset-1)] + SK.chemical_potential = ar['DMFT_results']['Iterations']['chemical_potential'+str(iteration_offset-1)].real + +iteration_offset = mpi.bcast(iteration_offset) +S.Sigma_iw = mpi.bcast(S.Sigma_iw) +SK.dc_imp = mpi.bcast(SK.dc_imp) +SK.dc_energ = mpi.bcast(SK.dc_energ) +SK.chemical_potential = mpi.bcast(SK.chemical_potential) + +# Calc the first G0 +SK.put_Sigma(Sigma_imp = [S.Sigma_iw]) +SK.calc_mu(precision=0.01) +S.G_iw << SK.extract_G_loc()[0] +SK.symm_deg_gf(S.G_iw, orb=0) + +#Init the DC term and the self-energy if no previous iteration was found +if iteration_offset == 0: + dm = S.G_iw.density() + SK.calc_dc(dm, U_interact=U, J_hund=J, orb=0, use_dc_formula=DC_type) + S.Sigma_iw << SK.dc_imp[0]['up'][0,0] + +mpi.report('%s DMFT cycles requested. Starting with iteration %s.'%(n_iterations,iteration_offset)) + +# The infamous DMFT self consistency cycle +for it in range(iteration_offset, iteration_offset + n_iterations): + + mpi.report('Doing iteration: %s'%it) + + # Get G0 + S.G0_iw << inverse(S.Sigma_iw + inverse(S.G_iw)) + # Solve the impurity problem + S.solve(h_int = H, **p) + if mpi.is_master_node(): + ar['DMFT_results']['Iterations']['Gimp_it'+str(it)] = S.G_iw + ar['DMFT_results']['Iterations']['Gtau_it'+str(it)] = S.G_tau + ar['DMFT_results']['Iterations']['Sigma_uns_it'+str(it)] = S.Sigma_iw + # Calculate double counting + dm = S.G_iw.density() + SK.calc_dc(dm, U_interact=U, J_hund=J, orb=0, use_dc_formula=DC_type) + # Get new G + SK.symm_deg_gf(S.Sigma_iw,orb=0) + SK.put_Sigma(Sigma_imp=[S.Sigma_iw]) + SK.calc_mu(precision=0.01) + S.G_iw << SK.extract_G_loc()[0] + + # print densities + for sig,gf in S.G_iw: + mpi.report("Orbital %s density: %.6f"%(sig,dm[sig][0,0])) + mpi.report('Total charge of Gloc : %.6f'%S.G_iw.total_density()) + + if mpi.is_master_node(): + ar['DMFT_results']['iteration_count'] = it + ar['DMFT_results']['Iterations']['Sigma_it'+str(it)] = S.Sigma_iw + ar['DMFT_results']['Iterations']['Gloc_it'+str(it)] = S.G_iw + ar['DMFT_results']['Iterations']['G0loc_it'+str(it)] = S.G0_iw + ar['DMFT_results']['Iterations']['dc_imp'+str(it)] = SK.dc_imp + ar['DMFT_results']['Iterations']['dc_energ'+str(it)] = SK.dc_energ + ar['DMFT_results']['Iterations']['chemical_potential'+str(it)] = SK.chemical_potential + +if mpi.is_master_node(): del ar diff --git a/doc/tutorials/images_scripts/nio.py.rst b/doc/tutorials/images_scripts/nio.py.rst new file mode 100644 index 00000000..489447c4 --- /dev/null +++ b/doc/tutorials/images_scripts/nio.py.rst @@ -0,0 +1,9 @@ +.. _nio.py: + +nio.py +----------- + +Download :download:`nio.py <./nio.py>`. + +.. literalinclude:: nio.py + :language: python diff --git a/doc/tutorials/images_scripts/plo.cfg b/doc/tutorials/images_scripts/plo.cfg new file mode 100644 index 00000000..91c72333 --- /dev/null +++ b/doc/tutorials/images_scripts/plo.cfg @@ -0,0 +1,27 @@ +[General] +BASENAME = nio +DOSMESH = -21 55 400 + +[Group 1] +SHELLS = 1 2 +EWINDOW = -9 2 +NORMION = False +NORMALIZE = True + +[Shell 1] # Ni d shell +LSHELL = 2 +IONS = 1 +CORR = True +TRANSFORM = 1.0 0.0 0.0 0.0 0.0 + 0.0 1.0 0.0 0.0 0.0 + 0.0 0.0 1.0 0.0 0.0 + 0.0 0.0 0.0 1.0 0.0 + 0.0 0.0 0.0 0.0 1.0 + +[Shell 2] # O p shell +LSHELL = 1 +IONS = 2 +CORR = False +TRANSFORM = 1.0 0.0 0.0 + 0.0 1.0 0.0 + 0.0 0.0 1.0 \ No newline at end of file diff --git a/doc/tutorials/images_scripts/plo.cfg.rst b/doc/tutorials/images_scripts/plo.cfg.rst new file mode 100644 index 00000000..1db06dae --- /dev/null +++ b/doc/tutorials/images_scripts/plo.cfg.rst @@ -0,0 +1,9 @@ +.. _plo.cfg: + +plo.cfg +----------- + +Download :download:`plo.cfg<./plo.cfg>`. + +.. literalinclude:: plo.cfg + :language: bash diff --git a/doc/tutorials/nio_csc.rst b/doc/tutorials/nio_csc.rst new file mode 100644 index 00000000..a2d7e94d --- /dev/null +++ b/doc/tutorials/nio_csc.rst @@ -0,0 +1,59 @@ +.. _nio_csc: + +DFT and projections +================================================== + +We will perform charge self-consitent DFT+DMFT calcluations for the charge-transfer insulator NiO. We still start from scratch and provide all necessary input files to do the calculations. + +VASP setup +------------------------------- +We start by running a simple VASP calculation to converge the charge density initially. +Use the :ref:`INCAR`, :ref:`POSCAR`, and :ref:`KPOINTS` provided and use your +own :file:`POTCAR` file. + +Let us take a look in the :file:`INCAR`, where we have to specify local orbitals as basis +for our many-body calculation. + +.. literalinclude:: images_scripts/INCAR + +`LORBIT = 14` takes care of optimizing the projectors in the energy window defined +by `EMIN` and `EMAX`. We switch off all symmetries with `ISYM=-1` since symmetries +are not implemented in the later DMFT scripts. Finally, we select the relevant orbitals +for atom 1 (Ni, d-orbitals) and 2 (O, p-orbitals) by the two `LOCPROJ` lines. +For details refer to the VASP wiki on the `LOCPROJ `_ flag. The projectors are stored in the file `LOCPROJ`. + + +plovasp +------------------------------ +Next, we postprocess the projectors, which VASP stored in the file `LOCPROJ`. +We do this by invoking :program:`plovasp plo.cfg` which is configured by an input file, e.g., named :ref:`plo.cfg`. + +.. literalinclude:: images_scripts/plo.cfg + +Here, in `[General]' we set the basename and the grid for calculating the density of +states. In `[Group 1]` we define a group of two shells which are orthonormalized with +respect to states in an energy window from `-9` to `2` for all ions simultanously +(`NORMION = False`). We define the two shells, which correspond to the Ni d states +and the O p states. Only the Ni shell is treated as correlated (`CORR = True`), i.e., +is supplemented with a Coulomb interaction later in the DMFT calculation. + +Converting to hdf5 file +------------------------------- +We gather the output generated by :program:`plovasp` into a hdf5 archive which :program:`dft_tools` is able to read. We do this by running :program:`python converter.py` on the script :ref:`converter.py`: + +.. literalinclude:: images_scripts/converter.py + +Now we are all set to perform a dmft calculation. + +DMFT +================================================== + +dmft script +------------------------------ + +Since the python script for performing the dmft loop pretty much resembles that presented in the tutorial on :ref:`srvo3`, we will not go into detail here but simply provide the script :ref:`nio.py`. Following Kunes et al. in `PRB 75 165115 (2007) `_ we use :math:`U=8` and :math:`J=1`. Here, we use :math:`\beta=5` instead of :math:`\beta=10` to speed up the calculations. + +Local lattice Green's function for all projected orbitals +---------------------- +We calculate the local lattice Green's function - now also for the uncorrelated orbitals, i.e., the O p states. Therefor we use :download:`NiO_local_lattice_GF.py AaAVFjc~u`WJXJ=;JJ$5-{ z?)S}oeDC{x-*mk`A~tYo$RoyfcTc)Ox+gVC(0{cOd)BRz8O%UI+sE-7PJhY1u`uvIop5 zJ}>q@J2SYa-?0-R0k@P@OY`0wLtN#Ka8;>Q4XfGdL41cY4`&Ay;`|7CkQ^QyiIa+g zovK2zSRRg7Up#crSoI2-xSx>~+@m=6)Kv_V@KI}o~x{uf!BKM9)7Zn zE|*Xce^*&mh-@x~Y_%d^#{U~4 z`a%$tcsM7XM>sI+GUUSd&f39GjhP3uE{KVN3E{!9B(zah(L@+w@69S+G2TVe6f>S; zT5nd-I66qzC}P4D9=x9kdd37(&>*F`6K!C9wGC|Y*ua)6IT-Bd98sX`p@LtKas8y5 zYJdsG3;$cm{XkgB48_f@O_ftB_yc$`$~7-#yD27 zW061xUXZPVAZPwzhcCRyVG=qnzocAeqa1Xpp+-D7bthMT>{W~ML7ON=CLfYrB-AIX zVgh_(3@{?a-k>9rrpPzxkg_7}J)~oUXOW7g0x!7hPIbf-c3$mN_ROnda-k|p1reVS zdT{sf1zZYFQwpTv8x%tdS#nZ#QF}rfUU9%hg)9aSP;mogICj7t&PH*=eo!-+aa68W zsXq@Dd|(!HS1t@^8;JL1v|!IWKQCKW!re05GJ<&V91r>#o^xP)Jem_f26w0$f+T$V z~`^Kd9==R7CE{8%f(A^{_ZsEe?Nh}y4yqX{iMltqBv37`I^+UbfXQ-My332%nWre$ zYsZ3DY&p!04?+vxy|P^Fd5a!J9lP-mk1PBj{|^?Ad5);xjq+F6b z!6cbYv!}@?b;q$as{dWP89VzmG1ZGF_1~pO5>k;ij9%Z`+n3)u@(x?Mx7T8jx~E^I z)5);#*9Z;1BsJbZPp#*GVVFn>@bXW26Bc6FDU>S^DU?B5up9F0SA|Mj>sEfA^*ed z{FC$nGzrB0-SP=iMob2H(5yn(`1ATfjs4bDs2w zehE)MxmATQrY9Z(J~*%wcxgIO?V`SCBp^bHt*ql42m)mbyGOev?dq85gkH~~>yz}_#NJQ$j?3`vc9sq|~qzB~dPj|w*fz#C^tvuA3JdkcI*OyDG70$XTgj+wmyF~2qEo0;dG5UN(i8}ib`x@iwnu)oz$RX}cGjvX7wxQdZ48gjcB(fUai&P@5L1D1|b!(8_7_)z2 z$AUxQ{@Z>U&JM#3Z!1Jwi`-&uOMLoTv>Nhva;TkqQasy^x+Gm-ugdq6uUE;lpH6;R zMm>U&XFu_5GGC9L26TFWx zjg_*ux(ZSE-XQtqR`f6s+q>0pd$?sJMM#m~4->0VnuLdm^Uz*0NTDvW;VD#&667F# z^pMYWa1gCTe06QYX16DS zunP`HJ|0oHqRbADQ7kauM(eR{vvG8BGYpkeAyB~*cLKch&I6aSNJxTj8$LxsN03L- zmLu}Y5j1sX_8}%ueqe|RDg>f^Vi*b&V*OG6c6;a~64hFjDU=7JcqoHnE{L^cha%n$n>BsN zr<3eOVyom+#O}eVh-Ayf(Pp5zwH?)v;6>uEIdzRmQ;zyNB*d2t6&E1BM9p7khV-CehR< zGL|m?4J34C=yzH|KU_?vRCCypVM9eph54j$K{L31odMJG%e|SMY8D$*Wb2J7Y{`}l zCcN?X;wF_e%5iQ7rTj0KWD|@Cn4ZRAdb>@dHM*vY&Ar}67Y34y z3#>l`B?6&1nutnm8}uQq7Mdga5Hdx{{Z|+^G1IAd)@0-)S9QQzkpjOQT~{34U`KDf zxwtff)pjmRNhA6S1lt{orjAq@Jf`e5(d-q2!w`#JO2^l0X8v)#mLs)iQGT?J|yl z&?x2BQtn9Cc~J1Bem^+Mh7foc4c5GZ`EW?q~T4qub384sE51-{MOWoA?f0xY3*U}f6N#Y z$Eaq&7-z?$ugZ!7vgb^9S(`CHIk@?ZH!V&B@Z7iTaDBoKHyt)$dX+?bnGqPJD((LG z$k`(@lxAd0uQ8SOeKfM958|5VJ*pBK*D`&&a{Px>v0op#cdO)JA80oDyeb~Ct1DeV zIcF$S9GxwW&J{-s`z)qFgI%GC=T%j*-Pg*yzESR^`igu-b5w@R8lT#4FWwhYzZJ7R z1;s>TN{^#@uj?q2ZY}T3T6fjE7R_%mi^z+=?=tZim|R~y4$f9be8iRqbbAm4>vFPi z1G+{TEYMyb#KRM8JxycMhwa-#g&{cyWqehC(@i-ZCJc8kY=aodQnFy4r4th#JXf;Q zAKnR9A_9)~S z##yl|#hWlvbhU2DG&6~8(Ys}ES&Ka%g$MUT)#QZh z46R{pMKP1+-7I?7W(-fE=9gT`&P!$(ttOL52@lTD!+rEio9av#cUyC8daY__nVUhh z^pN>yzy&jB80L^lIw_p{N{6Y}eE}nSn{xN^ST{3+q_;RaYN>|krwz*3R|aDY-uJZu z1EVC}i%Rofoo^h~S|7@l7+G%Xn zBp6tfnYHd_Dl=CcEp(XG5Lf1aWGOyDd}m|>rVJMrhk1sF+3z=UgX$Q|OR9nnA&x#c zo4lqv?xbG=+ zMg_ryZ*mJ_4WJN&L27l0MX*e^6B)3gA{C0Z{H+3Y5s%n22iyI7N!sIfwsXumIq!bo zz3;o<|IdpY8x(zSV9<@D1&zf&P>ETK^qD@a@L$^y@eGTHQt_O|SGt6qJH^73aK{O; zD3kA3wF|gP!d)k{#-lGI@VRv2Xh3>)SZ*k_2udekZN#^T3;a;9%4Jj>l2Fq@+ZKcn z{3@$RoAuPW`-y3{1;IPnX?Tmp4Su*!p9(#?6Iyf!By}BW@TKmOkOP93^=ryE3bjzu z0j;{;R52h&?T?iTOVEL#&V1~{P4=tiv5UUWEyA;&BHK2hZumC0fUY5@A+@|T?i>m* zJZA>XykU<2$JEp_J!a`=P+1tCJW*9fzu6f-6*YG7Z9IU^VA7INjhV$j8 z53aYTrGw_%A2xbc3g-N9W4ozKiTv=h+wMbcmVKBPboA!A4%Uq*!u2;X_B;hTs(|q$ zT_n6#GmwYj0tGs*yn(G4)5-|SH0juexa5p$MlzW}A3nHVFdOei2K+PEl!vg+A>kUq zo8ra3TZNXbz_bF_GLjXjVHBI+guXzi3+9gAB0a8}?B1?j0gk->maj!nn0)!B#?_Ca zWc*un$>gag(_Xf58>qoIsoAdP5=F)^iOk&wx`kPx>fZ7cLB{WjaCZt8@!`9utut;G z94AIoh5KOd(>}msM5NGwswdHp(E6smx9N)kPaTY3K(M1~)cAMBJC;jXQu?LuSyGlr z*h)tZ30Ji@9@&K;IHI4XDs>pfo*WLUfw@It20mdG%bgHtai|`rLaw*rhYZi31*d)zdIq`mJ(BOH z+{IZK0`xGK!_p*P?8dM91L|{ zApQ$9@~#hxEKY(aQ3SVs(Aas4U^snNGVHV%=CvttM3APHimE#BiGZw&EzkTORPAiz z+UPh$Jrf3eBpq|1+IHMhGpZuvPid}DGRZV;PZtAN^t7MMUycTW^YBFa_5+7y`AjZH z)vSDiuYGv`eu`ya_fhI@!>Xk_1f|6+kiAceUL*DA|VrnL>*k+(##>QW3WBy$41kZ`MJt)s4i&+9p zw!3i|h9SoqSn`?&V^2+|I$wUOs2`J+3i|8Esg34PUF5-}(<|xTvC~Dp^#JiF%FzHT z(U1B?H~DI*X8C8c7J-f~9Ri=vgW^?#AaTM0`k6&+ZPc|2z5AQiotYY(-r6ZR>+6mc zlHvap`IbOdn({cMgw8b?GuUTUF1Y%|R0%Fk+--2~|E*O5$0m6f(XJ{bXGnllrA-KG;hg zBySmdmL~cam!aMk7zf8&29iO~qeARi5*zzGT1vY?@~7oU1Wt?wPE-XmM_*|;SZ_(9 zA;z+IL+iOLE^*yzA@O-cii2oy)_%bP)eU9rCFDa6lDifSZcKEc%dN8I6=#%^#i9w0 z3=k3MCwDcV8hZc!26;cw9II-IPR4{$5i-eK_HdHt!st%yt|8PrH?o|p28FUkwPlGd z*0TBjcSp>7*SDKB4sMtkf_d<6gZ{2o18eSS-q6% zPH$Zr*Np`v%9)wBhBVip;c(X+2Qd@K-xyCK){T1jBDfE^)!hj&$UUohupXyZ6fAHWY|9!M<{LBar;zgUOtU{`++4Y?JPn`cdY8K zXiqMgdl*B@t(HKnTeT+nP+eEd9uZ$SNGFDubH5LK^M!x#7o_q?0W%C1kc@hWRi~Okw*{9+ekPY2h@5UH+GIjBxuJ&q4w!}#i80w znTL=sR?r`pBl2o1dK{|Sb-PpcrhhdS&($_Qd5xx6rK&{EY?dwB()n@z&QZgel%BoVhqKjoo+P|y5zi}k(hE?xry4(r=>k9A-TGQQAufVE%g$9pI*iP~` zYpKp6qB>sWdDyp#_Cl-p>$D+N+kzZR6NiC+$FYoPSBNkDBMs8HNptaDk`%tKS&hHi z16oEh%1%!@%x8W~$5W?TO^-a{dLWc^ig3?^HE&PDYg6yB$k5Xu&DC-<>pe~4 z8kT>__7KX2>f6fktGzgNXM){*HcQf*qgHg2pRX4I{@BV5_VkR#6ho@cBGYJnEmJJc zP<*U#20AVc7)*RywGy?G5~WsD%eHbxgjQ&#QOjQEa?!y~t=h_$-qi-)#8*u2lFfV> zq{s+gg>N+|!(8?|%dsFQ=ODGBH|XyL!os1}e3zXumBLj|&w#c`_k-O#hE9sjeY zlxFT?kar6=bA!J(O_}S%?`V{x;0TB6NHBjQJcYc^U5llV${ifbQV4Y3;$jmmoxy*7 zYmr*f%)O2k)OZH#^Q>zuCrqR_SypIGgS6AhPO0u~mRue(dOi|j^1-3lf8hTU8qrd> zPop2h$gzEzXi8{P5J7FEaX(EGT_0L3aq8KZ6qN7PgmN`_s=arw={PJ9!&Q1A7)S?at$I!3~^CAN3 zCtPaUjltzdECpb9YH*;IB6(d_^`y;>;p~=f*4UVZNtk|3!n8kP`8k}?DsoQ>lE!n? zlfvh@r&QfzkZIl2k63CBBss3-9s~cKt%(0En(WRoKVskSaJt23v6zPF$2bBPL{*_* z<1%)bPbD%rplu^8jSSXh1IUBK2z-;lJ@{R^>b#&qX=PMW;&vt-8@0`3?l+Wo{UYa+ zeX&ZAMuzG72gk73ZX>9NVyZ51wv+QS*=Fo?r1_@dKz@j9mWI`2dW#1`8f5ODYDT{n zf!?7DQwtf&uz+-C0VV_Lcrj$#B3o`{CGW9?!?Jj8A-ULAq_$?+SUw*kROspZvfL3&;=Ss>uE=xJflwigdlOv<~W?@VmhHW72P#-3L aG;H^)H|Jj!bJ=ab$->dx*v4F&-SKbnBt3cn diff --git a/test/plovasp/converter/pg_output.out.h5 b/test/plovasp/converter/pg_output.out.h5 index bb0da56fd93ad43d15fb53e8fe130c462f643a13..5427ce6cb70db86aa24db40328d61537a333c110 100644 GIT binary patch delta 5754 zcma(VZERE5^}Tl!^6ZfIIU~dgU^^ux21@3900sqcpaCZZ;sB)q3U92WInjG?jNN7{KA6w8?Dg84={V~NA9cf4DtLV0D9n$W7=R7~Z zM?7#8-*?VC_k5po&-Kjj>-6_MxgbdmR;e4WPcHiAMO9a5#HY|H>4vyWi|-=ML(PKm ziNYc*G*lT+!zR279s&p#67W!vU>#op^oj^@Sef|rHzE@vk7TL=wTPsltw1UYRRX8> z!fQccfBd7j0U`wpqgk3RW75%2s`u?lodt4~FHC(2;^rQ}8bofP7l`p>(hVb=ZY65* zW2(2lKnH3%Yd8%N!}L2(OU>XP)|?cCtXia#Z^*Z*>$v)b9yDUI;ki4^pBb7^meth= zGsssz$Jni@O(B^7LSo?6Ids8%qoIRgK7X3o0kDj|m;q1_df>rScDxth={>N{BEKxU z*rW@64liQ&VV!}62sy7>^tlSievE$7<)_NvWi&U!&8_0*(m?$^L`IiF#+E`3TgX_$ z15uz?YT7L4T{nHpUi`yT{-UK)ghOa`~l{tAN1kn<)B3=7;@gM#j=uS(}#_(jlA)7OwG*GxOG zu$&n2d8NWeI?$f*FWms9w!)Sao$$4}0v%2UJPN6H9)%|PC3smX5QJlN$n%=2&!7(v ztV5@%rFCMp&QlJ*LBg734tuFX_OEpBS+uqAt_@0|6KB z1IvXli}d3Y=cW9%iQ*Qqzf!5@_JXjdHEsPBq@y2r4x_wqrg z{Z$}=N}L}8$TE>@GIQ|%lyb(#Abeyg2H7Dx(0j2tFbFTrNJbG(4ozTlML*{3)=I&{ zJcG5+R^nPw+5@~w7FTR&R=zlPJhP9|CtznvyXGT7`8Yrh9}UwlkAH^xkGYp9!@x0a zgpq5`ld5C{@3GpSy2&OQIJN(IE{BN0TPiqLuH%O-z8@14CBD-^-sgMY=T6Bjo&Y>SL z$YqY389EUSsQMpSx2>n!&mBT5g}NuLOsmB!hLv90d;UESCm^EoFCZ-37wPqyMq95# zbJ?kiE}1LqPpnf`IE`d%4k>M99|?$O%Ia?3#$Zg-|0XqwXX$|Vl$fI8UcbjVtO`H9 zJNmP`qvGwmqf*wv%gK}hsM>-m%SWihqXQQ-fsBP>gXI(mmN-0Z7Hez2l z(k$`Nk!Jk1;tqi{Vex_X?jtE7y0k}+B-O>SM9hmsyZ)`@_wi%Pyh`N120l2;F+z4k zZ{ef*W=RO}jqH<77)xmwANve8P2EE!;I$C@N}4FviPh$r2zjDS|5N&*!bylIJO`u& zk^Z8_TV7@SY;f|mIYljQ{#A5{D5}2f_NjkMH&^qZ6G4}7rd(&hX`QLC8Onbc!^vn2 zsWx`>yS{ZRt(%QOaO6-QJ9{Mco#w`)c#4`mbREF_I&jN>*SNOPRJB~2MuNirD(&yR zR0`IfGSdU~V=OzpB)gvQQrWlkteWlzF)5T&=_EDHPSFqCjdXZ#i(C(WD631BOHM28 z!8hiNP7rAAi?w)~9ZsEn35YU1WTLpnxS61d(d}g>@;TJ$|B}S3&pIS1y>=Ho>%QSW z`##f!xL^Ieg-&-VWO@tFA3m>>*{%GM4eTV>WgywaP70fiM|OcC97++p+s6$iJ`oAOFJF&o$;LPd`&hhrqEAIwqI-e2vx-Vt!FjWu z#w$ad(bQfNZ+YwTYnMih;2Gzt!I^($j)N37j(nfa$}&|Hh2}kG4Q2<*TwEB;lk9Gm ztgH$0t9R~@LqjdoNiL&7xw4&Fpo^qb;ylQu46ei^(fKB3&QBso~-(_LE z-8V#tQW~?`Xo1Bs77rL}CjK(oqF*BP97g?s7N!HmQDU;^xQiY$IB?T zbF&##>Zgl6pnr;flx(VfQe%F~4hO-XJh{Uws*L1*yW6e(=UZv{a@Kv(XA39z-4Lfd zA?WyN&UuMSV}6e}@hR5)S#udFdaQZv%Meo%ykm;PlFFhjl!ag@O%;$g<|m!63g_e ztR!ESe8#7DQi2A_+yzra+mXSL`6={K;SrZK9;G7he z!2}H9R!IIOwxi-*Ib~M!5QS z@HkqFEC#rF)-M1u3GZ;bvHxy(!zNdDJy$s(YEX4P! zd~bSvw)hs3`Z(-z%*luYVqBDK;yJ-n_3>w>-( zP6km)z*rd`u8HJ1o2W6rla>O>DBO93;p2VW@QDUlK5o~T`)j%b<>wwrO*lr4nvTi$ z0Vc)|sUTd9SdKx{%-aI9xSp3|G2sB=p3wozJqn9P!9C$^=pI>nL9DoYqyqavG?*>N zYmBhq(B%1xjH>lE7hRj&>!KPEKU*8n?Q(#-O3kkDm!lDA`oI3Ns}1Vz@)uH7YpBlm z|Ez%*5?BZmjC-LcKUB{PP9L$_5pQLGg51jNdW4IKi+J+rQi;ch%*4%fshtbE;Iid) z^(BBKLq$wSRJ+vdm{=IoL4*sosqEl)%9IBq)Tlk?P6B3v-3q>HnR}JPH?B`-%PP~U z7WMYZMkFH%izN5MK;-~H?%X0IhU!~>P!V>I1%*zf*}@$CHjifiB;w^+%NgD5q0%Y7 zAI;%@^%dZUAkhv392TkpOA))~_+XW2`oI!aYDo5jR&ua}NXAfKumF3|McJt5~_IdCC%9pPR ziNeT7(#t|Jmf7oX!xmXS^tgg47{ee~`L(48CR)+>82fzqI?7dW$(5H2W2Ru~Zl&Gp z0NZ4(uO51NL?RNgj8WZfUgE8=XMAF*;#7}|>=mxGd(4evA-5-`QKQ)A3l_^7?|NJ` z@2F&pcDX8tnL06oq~h$v1WrEbiFcZobD0TInRN1oU#zH9xj>Ao=V{Qz?nOH-l>khW zqYq&qu=%jasZi{V$5XWoXvCc$!V{IVW^Lhi8<^VwoGg=CL%27->Jt;KDfxzIO&wv+ z033cD*pJTK`l}o?hX!>0{L@ag5LL>?M7==S+7_Kj0IkAV_*&-JgHl3VG!d&LsHb^? zuP9+g$M*9PeBJxg2r769K{D>i&{^O?QkWnXUe1#+8H+rc&cgG~6H{;?lM-ub>6Wz@ zv&Fr$ zAvJZlaoa0@2|;D%@jz+;Aa4Lz)23B!Q)m`fH1#wp>5En5fwvMrymIbGSI$W@SI)`Z z-)Gftb;!46_UGY`P;&l`v4?}gNM26fo{4Vv@k%AK&O3oy?jIIo?LgtH`U*T=1oL#D ziMukp(R;eVl#oh-4vB0rI1}6EOG-u|=KE&bg^Q5w!611LGfeWl_&m@+U5`P#P1M90e!&nwqk zf2`oBJ`iqMWsW~I`A~T90{i!gfDx|BBWCtdNm~U(3ljTDlLqr1f03AWDL(|aHf(+= zEJxhfn`V9~(d%nXA6NeM1ZHZ;OmjTDS_~|!YaZ>U4r7P+>>)&EKRwbJTMKp4N{i|~ zfMaMuD^1D>yV8We#FSkHJekICbqFXlN} zh*Gy=3c+-*05@1N=XP5h;d~ZhrJtjHCMJDkP;nTT@123jdbnxD!P2^f$gx?S)};HK z!+`O0r+SDMi_b&0PjvyWybfPGVOpTfFDGGCFK+y;)n~psX`)AHDxRHG-()WNyOR!f zn@j%iq>;_$QoKBsbNKQ^Z5L{U;SU?!t=8RLth>uwyIs8;3QFBwweGIw+HLe;$w$E? z38Q#N4#y&Hrz|xjHM-bW>lD342H}+WHyo#wo&{cWeIBQX9~o%4U(mT=93P`+(R^Wt zn33B@Cq??xbiYaQ@;>@DqDAVV&22@-uld{d0KX)rn;i_xzxVdhi)SX;do72M6dLxA zNBHk+juX(ZxZB9Em0iFFTH@4Bfc~#&zr3Gz$UExT)f7S4HMPUB>-+$)B0Jz~>({Vo zEsoTRS{Qx;QHAfMY2)S)J)JEpM(Ln0!++QL`O~W#GGA~p&ph`e%p*hoHm3lYDQ#Vs`k!y0zs=?k9vYDiFGBlPSlDO2?X^u0(?J(VCgRO zB8!T5IGsd1?Q1ulx2{>8BZH}B(etJnlFNRwpZ|&+HY#7r{r=c1poz+U+|)sz=~}~H d3S?*0LsFc=fAPyrb*0Zgiipz8z2lUs{{?*eI^qBT