11from itertools import *
2- import numpy as np
32import triqs .utility .mpi as mpi
4- from h5 import *
3+ from h5 import HDFArchive
54from triqs .gf import *
65import sys , triqs .version as triqs_version
76from triqs_dft_tools .sumk_dft import *
@@ -25,13 +24,11 @@ def dmft_cycle():
2524 Converter .convert_dft_input ()
2625
2726 beta = 5.0
28- mesh = MeshImFreq (beta = beta , S = 'Fermion' , n_iw = 1000 )
27+ mesh = MeshImFreq (beta = beta , S = 'Fermion' , n_iw = 500 )
2928
3029 SK = SumkDFT (hdf_file = filename + '.h5' , use_dft_blocks = False , mesh = mesh )
3130
32- Sigma = SK .block_structure .create_gf (mesh = mesh )
33- SK .put_Sigma ([Sigma ])
34- G = SK .extract_G_loc ()
31+ G = SK .extract_G_loc (transform_to_solver_blocks = False , with_Sigma = False )
3532 SK .analyse_block_structure_from_gf (G , threshold = 1e-2 )
3633 for i_sh in range (len (SK .deg_shells )):
3734 num_block_deg_orbs = len (SK .deg_shells [i_sh ])
@@ -52,7 +49,7 @@ def dmft_cycle():
5249 mpi .report ('GF struct sumk: %s' % SK .gf_struct_sumk )
5350 mpi .report ('GF struct solver: %s' % SK .gf_struct_solver )
5451
55- S = Solver (beta = beta , gf_struct = gf_struct , n_iw = 1000 )
52+ S = Solver (beta = beta , gf_struct = gf_struct , n_iw = 500 )
5653
5754 # Construct the Hamiltonian and save it in Hamiltonian_store.txt
5855 H = Operator ()
@@ -75,14 +72,16 @@ def dmft_cycle():
7572 p = {}
7673 p ["max_time" ] = - 1
7774 p ["random_name" ] = ""
78- p ["random_seed" ] = 123 * mpi .rank + 567
79- p ["length_cycle" ] = 100
75+ p ["length_cycle" ] = 400
8076 p ["n_warmup_cycles" ] = 2000
81- p ["n_cycles" ] = 20000
77+ p ["n_cycles" ] = 80000
8278 p ["fit_max_moment" ] = 4
83- p ["fit_min_n " ] = 30
84- p ["fit_max_n " ] = 50
79+ p ["fit_min_w " ] = 20
80+ p ["fit_max_w " ] = 30
8581 p ["perform_tail_fit" ] = True
82+ p ["measure_density_matrix" ] = True
83+ p ["use_norm_as_weight" ] = True
84+
8685
8786 # Double Counting: 0 FLL, 1 Held, 2 AMF
8887 DC_type = 0
@@ -93,32 +92,32 @@ def dmft_cycle():
9392
9493 iteration_offset = 0
9594 if mpi .is_master_node ():
96- ar = HDFArchive (filename + '.h5' , 'a' )
97- if not 'DMFT_results' in ar :
98- ar .create_group ('DMFT_results' )
99- if not 'Iterations' in ar ['DMFT_results' ]:
100- ar ['DMFT_results' ].create_group ('Iterations' )
101- if not 'DMFT_input' in ar :
102- ar .create_group ('DMFT_input' )
103- if not 'Iterations' in ar ['DMFT_input' ]:
104- ar ['DMFT_input' ].create_group ('Iterations' )
105- if not 'code_versions' in ar ['DMFT_input' ]:
106- ar ['DMFT_input' ].create_group ('code_versions' )
107- ar ['DMFT_input' ]['code_versions' ]["triqs_version" ] = triqs_version .version
108- ar ['DMFT_input' ]['code_versions' ]["triqs_git" ] = triqs_version .git_hash
109- ar ['DMFT_input' ]['code_versions' ]["cthyb_version" ] = cthyb_version .version
110- ar ['DMFT_input' ]['code_versions' ]["cthyb_git" ] = cthyb_version .triqs_cthyb_hash
111- ar ['DMFT_input' ]['code_versions' ]["dft_tools_version" ] = dft_tools_version .version
112- ar ['DMFT_input' ]['code_versions' ]["dft_tools_git" ] = dft_tools_version .triqs_dft_tools_hash
113- ar ['DMFT_input' ]['sumk_block_structure' ] = SK .block_structure
114- if 'iteration_count' in ar ['DMFT_results' ]:
115- iteration_offset = ar ['DMFT_results' ]['iteration_count' ]+ 1
116- S .Sigma_iw = ar ['DMFT_results' ]['Iterations' ]['Sigma_it' + str (iteration_offset - 1 )]
117- SK .dc_imp = ar ['DMFT_results' ]['Iterations' ]['dc_imp' + str (iteration_offset - 1 )]
118- SK .dc_energ = ar ['DMFT_results' ]['Iterations' ]['dc_energ' + str (iteration_offset - 1 )]
119- SK .chemical_potential = ar ['DMFT_results' ]['Iterations' ]['chemical_potential' +
120- str (iteration_offset - 1 )].real
121- ar ['DMFT_input' ]["dmft_script_it" + str (iteration_offset )] = open (sys .argv [0 ]).read ()
95+ with HDFArchive (filename + '.h5' , 'a' ) as ar :
96+ if 'DMFT_results' not in ar :
97+ ar .create_group ('DMFT_results' )
98+ if 'Iterations' not in ar ['DMFT_results' ]:
99+ ar ['DMFT_results' ].create_group ('Iterations' )
100+ if 'DMFT_input' not in ar :
101+ ar .create_group ('DMFT_input' )
102+ if 'Iterations' not in ar ['DMFT_input' ]:
103+ ar ['DMFT_input' ].create_group ('Iterations' )
104+ if not 'code_versions' not in ar ['DMFT_input' ]:
105+ ar ['DMFT_input' ].create_group ('code_versions' )
106+ ar ['DMFT_input' ]['code_versions' ]["triqs_version" ] = triqs_version .version
107+ ar ['DMFT_input' ]['code_versions' ]["triqs_git" ] = triqs_version .git_hash
108+ ar ['DMFT_input' ]['code_versions' ]["cthyb_version" ] = cthyb_version .version
109+ ar ['DMFT_input' ]['code_versions' ]["cthyb_git" ] = cthyb_version .triqs_cthyb_hash
110+ ar ['DMFT_input' ]['code_versions' ]["dft_tools_version" ] = dft_tools_version .version
111+ ar ['DMFT_input' ]['code_versions' ]["dft_tools_git" ] = dft_tools_version .triqs_dft_tools_hash
112+ ar ['DMFT_input' ]['sumk_block_structure' ] = SK .block_structure
113+ if 'iteration_count' in ar ['DMFT_results' ]:
114+ iteration_offset = ar ['DMFT_results' ]['iteration_count' ]+ 1
115+ S .Sigma_iw = ar ['DMFT_results' ]['Iterations' ]['Sigma_it' + str (iteration_offset - 1 )]
116+ SK .dc_imp = ar ['DMFT_results' ]['Iterations' ]['dc_imp' + str (iteration_offset - 1 )]
117+ SK .dc_energ = ar ['DMFT_results' ]['Iterations' ]['dc_energ' + str (iteration_offset - 1 )]
118+ SK .chemical_potential = ar ['DMFT_results' ]['Iterations' ]['chemical_potential' +
119+ str (iteration_offset - 1 )].real
120+ ar ['DMFT_input' ]["dmft_script_it" + str (iteration_offset )] = open (sys .argv [0 ]).read ()
122121 iteration_offset = mpi .bcast (iteration_offset )
123122 S .Sigma_iw = mpi .bcast (S .Sigma_iw )
124123 SK .dc_imp = mpi .bcast (SK .dc_imp )
@@ -151,10 +150,11 @@ def dmft_cycle():
151150 # Solve the impurity problem
152151 S .solve (h_int = H , ** p )
153152 if mpi .is_master_node ():
154- ar ['DMFT_input' ]['Iterations' ]['solver_dict_it' + str (it )] = p
155- ar ['DMFT_results' ]['Iterations' ]['Gimp_it' + str (it )] = S .G_iw
156- ar ['DMFT_results' ]['Iterations' ]['Gtau_it' + str (it )] = S .G_tau
157- ar ['DMFT_results' ]['Iterations' ]['Sigma_uns_it' + str (it )] = S .Sigma_iw
153+ with HDFArchive (filename + '.h5' , 'a' ) as ar :
154+ ar ['DMFT_input' ]['Iterations' ]['solver_dict_it' + str (it )] = p
155+ ar ['DMFT_results' ]['Iterations' ]['Gimp_it' + str (it )] = S .G_iw
156+ ar ['DMFT_results' ]['Iterations' ]['Gtau_it' + str (it )] = S .G_tau
157+ ar ['DMFT_results' ]['Iterations' ]['Sigma_uns_it' + str (it )] = S .Sigma_iw
158158 # Calculate double counting
159159 dm = S .G_iw .density ()
160160 SK .calc_dc (dm , U_interact = U , J_hund = J , orb = 0 ,
@@ -171,20 +171,21 @@ def dmft_cycle():
171171 mpi .report ('Total charge of Gloc : %.6f' % S .G_iw .total_density ().real )
172172
173173 if mpi .is_master_node ():
174- ar ['DMFT_results' ]['iteration_count' ] = it
175- ar ['DMFT_results' ]['Iterations' ]['Sigma_it' + str (it )] = S .Sigma_iw
176- ar ['DMFT_results' ]['Iterations' ]['Gloc_it' + str (it )] = S .G_iw
177- ar ['DMFT_results' ]['Iterations' ]['G0loc_it' + str (it )] = S .G0_iw
178- ar ['DMFT_results' ]['Iterations' ]['dc_imp' + str (it )] = SK .dc_imp
179- ar ['DMFT_results' ]['Iterations' ]['dc_energ' + str (it )] = SK .dc_energ
180- ar ['DMFT_results' ]['Iterations' ]['chemical_potential' + str (it )] = SK .chemical_potential
174+ with HDFArchive (filename + '.h5' , 'a' ) as ar :
175+ ar ['DMFT_results' ]['iteration_count' ] = it
176+ ar ['DMFT_results' ]['Iterations' ]['Sigma_it' + str (it )] = S .Sigma_iw
177+ ar ['DMFT_results' ]['Iterations' ]['Gloc_it' + str (it )] = S .G_iw
178+ ar ['DMFT_results' ]['Iterations' ]['G0loc_it' + str (it )] = S .G0_iw
179+ ar ['DMFT_results' ]['Iterations' ]['dc_imp' + str (it )] = SK .dc_imp
180+ ar ['DMFT_results' ]['Iterations' ]['dc_energ' + str (it )] = SK .dc_energ
181+ ar ['DMFT_results' ]['Iterations' ]['chemical_potential' + str (it )] = SK .chemical_potential
181182
182183 if mpi .is_master_node ():
183184 print ('calculating mu...' )
184185 SK .chemical_potential = SK .calc_mu (precision = 0.000001 )
185186
186187 if mpi .is_master_node ():
187- print ('calculating GAMMA ' )
188+ print ('calculating charge density update ' )
188189 SK .calc_density_correction (dm_type = 'vasp' )
189190
190191 if mpi .is_master_node ():
@@ -196,10 +197,8 @@ def dmft_cycle():
196197 SK .calc_dc (dm , U_interact = U , J_hund = J , orb = 0 , use_dc_formula = DC_type , use_dc_value = DC_value )
197198
198199 if mpi .is_master_node ():
199- ar ['DMFT_results' ]['Iterations' ]['corr_energy_it' + str (it )] = correnerg
200- ar ['DMFT_results' ]['Iterations' ]['dc_energy_it' + str (it )] = SK .dc_energ [0 ]
201-
202- if mpi .is_master_node ():
203- del ar
200+ with HDFArchive (filename + '.h5' , 'a' ) as ar :
201+ ar ['DMFT_results' ]['Iterations' ]['corr_energy_it' + str (it )] = correnerg
202+ ar ['DMFT_results' ]['Iterations' ]['dc_energy_it' + str (it )] = SK .dc_energ [0 ]
204203
205204 return correnerg , SK
0 commit comments