diff --git a/python/solid_dmft/dmft_cycle.py b/python/solid_dmft/dmft_cycle.py index 1ef26251..9a6fd41e 100755 --- a/python/solid_dmft/dmft_cycle.py +++ b/python/solid_dmft/dmft_cycle.py @@ -57,6 +57,7 @@ from solid_dmft.dmft_tools import manipulate_chemical_potential as manipulate_mu from solid_dmft.dmft_tools import initial_self_energies as initial_sigma from solid_dmft.dmft_tools import greens_functions_mixer as gf_mixer +from solid_dmft.io_tools.dict_to_h5 import prep_params_for_h5 def _extract_quantity_per_inequiv(param_name, n_inequiv_shells, general_params): """ @@ -370,7 +371,7 @@ def dmft_cycle(general_params, solver_params, advanced_params, dft_params, map_imp_solver.append(isolver) break solver_type_per_imp = [solver_params[map_imp_solver[iineq]]['type'] for iineq in range(sum_k.n_inequiv_shells)] - mpi.report('DEBUG', solver_type_per_imp) + mpi.report(f'Solver type per impurity: {solver_type_per_imp}') # Checks that enforce_off_diag true for ftps and hartree if any(s in ['ftps', 'hartree'] and not e for s, e in zip(solver_type_per_imp, general_params['enforce_off_diag'])): @@ -467,7 +468,7 @@ def dmft_cycle(general_params, solver_params, advanced_params, dft_params, archive['DMFT_input']['rot_mat'] = sum_k.rot_mat else: previous_rot_mat = None - if deg_orbs_ftps is not None: + if 'solver_struct_ftps' in archive['DMFT_input']: deg_orbs_ftps = archive['DMFT_input/solver_struct_ftps'] sum_k.block_structure = mpi.bcast(sum_k.block_structure) @@ -507,10 +508,10 @@ def dmft_cycle(general_params, solver_params, advanced_params, dft_params, # If new calculation, writes input parameters and sum_k <-> solver mapping to archive if iteration_offset == 0: if mpi.is_master_node(): - # FIXME: writing to archive crashes because of parameters that are None - # archive['DMFT_input']['general_params'] = general_params - # archive['DMFT_input']['solver_params'] = solver_params - # archive['DMFT_input']['advanced_params'] = advanced_params + archive['DMFT_input']['general_params'] = prep_params_for_h5(general_params) + archive['DMFT_input']['solver_params'] = prep_params_for_h5(solver_params) + archive['DMFT_input']['dft_params'] = prep_params_for_h5(dft_params) + archive['DMFT_input']['advanced_params'] = prep_params_for_h5(advanced_params) archive['DMFT_input']['block_structure'] = sum_k.block_structure archive['DMFT_input']['deg_shells'] = sum_k.deg_shells diff --git a/python/solid_dmft/dmft_tools/solver.py b/python/solid_dmft/dmft_tools/solver.py index c666fcc9..abbf8751 100755 --- a/python/solid_dmft/dmft_tools/solver.py +++ b/python/solid_dmft/dmft_tools/solver.py @@ -32,6 +32,8 @@ import triqs.utility.mpi as mpi from h5 import HDFArchive +from solid_dmft.io_tools.dict_to_h5 import prep_params_for_h5 + from . import legendre_filter from .matheval import MathExpr @@ -142,10 +144,12 @@ def __init__(self, general_params, solver_params, advanced_params, sum_k, icrsh, self.h_int = h_int self.iteration_offset = iteration_offset self.deg_orbs_ftps = deg_orbs_ftps - if solver_params.get("random_seed") is None: + + # Stores random_seed as MathExpr object to evaluate it at runtime + if self.solver_params.get('random_seed') is None: self.random_seed_generator = None else: - self.random_seed_generator = MathExpr(solver_params["random_seed"]) + self.random_seed_generator = MathExpr(self.solver_params['random_seed']) # initialize solver object, options are cthyb if self.solver_params['type'] == 'cthyb': @@ -261,11 +265,9 @@ def _init_ImFreq_objects(self): self.Delta_time = self.G_time.copy() # create all Legendre instances - if (self.solver_params['type'] == 'cthyb' and self.solver_params['measure_G_l'] - or self.solver_params['type'] == 'cthyb' and self.general_params['legendre_fit'] - or self.solver_params['type'] == 'ctseg' and self.solver_params['measure_gl'] - or self.solver_params['type'] == 'ctseg' and self.general_params['legendre_fit'] - or self.solver_params['type'] == 'hubbardI' and self.solver_params['measure_G_l']): + if (self.solver_params['type'] in ['cthyb', 'ctseg'] + and (self.solver_params['measure_G_l'] or self.solver_params['legendre_fit']) + or self.solver_params['type'] == 'hubbardI' and self.solver_params['measure_G_l']): self.n_l = self.solver_params['n_l'] self.G_l = self.sum_k.block_structure.create_gf(ish=self.icrsh, gf_function=Gf, space='solver', @@ -347,13 +349,12 @@ def solve(self, **kwargs): solve impurity problem with current solver ''' - if self.random_seed_generator is None: - random_seed = {} + if self.random_seed_generator is not None: + self.triqs_solver_params['random_seed'] = int(self.random_seed_generator(it=kwargs["it"], rank=mpi.rank)) else: - random_seed = {"random_seed": int(self.random_seed_generator(it=kwargs["it"], rank=mpi.rank))} + assert 'random_seed' not in self.triqs_solver_params if self.solver_params['type'] == 'cthyb': - if self.solver_params['delta_interface']: mpi.report('\n Using the delta interface for cthyb passing Delta(tau) and Hloc0 directly.') # prepare solver input @@ -370,7 +371,6 @@ def solve(self, **kwargs): else: self.triqs_solver.Delta_tau[name] << make_gf_from_fourier(self.Delta_freq[name], self.triqs_solver.Delta_tau.mesh, tail) - # Make non-interacting operator for Hloc0 Hloc_0 = Operator() for spin, spin_block in solver_eal.items(): @@ -393,13 +393,12 @@ def solve(self, **kwargs): if not 'it_-1' in archive['DMFT_input/solver']: archive['DMFT_input/solver'].create_group('it_-1') archive['DMFT_input/solver/it_-1'][f'S_{self.icrsh}'] = self.triqs_solver - # TODO: write only params that the solver can use - archive['DMFT_input/solver/it_-1'][f'triqs_solver_params_{self.icrsh}'] = self.triqs_solver_params + archive['DMFT_input/solver/it_-1'][f'triqs_solver_params_{self.icrsh}'] = prep_params_for_h5(self.triqs_solver_params) archive['DMFT_input/solver/it_-1']['mpi_size'] = mpi.size # Solve the impurity problem for icrsh shell # ************************************* - self.triqs_solver.solve(h_int=self.h_int, **{ **self.triqs_solver_params, **random_seed }) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # ************************************* # call postprocessing @@ -409,13 +408,13 @@ def solve(self, **kwargs): # fill G0_freq from sum_k to solver self.triqs_solver.G0_iw << self.G0_freq - if self.general_params['h_int_type'] == 'dynamic': + if self.general_params['h_int_type'][self.icrsh] == 'dynamic': for b1, b2 in product(self.sum_k.gf_struct_solver_dict[self.icrsh].keys(), repeat=2): self.triqs_solver.D0_iw[b1,b2] << self.U_iw[self.icrsh] # Solve the impurity problem for icrsh shell # ************************************* - self.triqs_solver.solve(h_int=self.h_int, **{ **self.solver_params, **random_seed }) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # ************************************* # call postprocessing @@ -425,14 +424,10 @@ def solve(self, **kwargs): # fill G0_freq from sum_k to solver self.triqs_solver.G0_iw << self.G0_freq - print() - # Solve the impurity problem for icrsh shell # ************************************* # this is done on every node due to very slow bcast of the AtomDiag object as of now - self.triqs_solver.solve(h_int=self.h_int, calc_gtau=self.solver_params['measure_G_tau'], - calc_gw=True, calc_gl=self.solver_params['measure_G_l'], - calc_dm=self.solver_params['measure_density_matrix']) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # if density matrix is measured, get this too. Needs to be done here, # because solver property 'dm' is not initialized/broadcastable if self.solver_params['measure_density_matrix']: @@ -450,9 +445,7 @@ def solve(self, **kwargs): # Solve the impurity problem for icrsh shell # ************************************* # this is done on every node due to very slow bcast of the AtomDiag object as of now - self.triqs_solver.solve(h_int=self.h_int, with_fock=self.solver_params['with_fock'], - one_shot=self.solver_params['one_shot'], - method=self.solver_params['method'], tol=self.solver_params['tol']) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # call postprocessing self._hartree_postprocessing() @@ -490,7 +483,7 @@ def make_positive_definite(G): self.Delta_freq_solver = make_positive_definite(self.Delta_freq_solver) # remove off-diagonal terms - if self.general_params['diag_delta']: + if self.solver_params['diag_delta']: for name, delta in self.Delta_freq_solver: for i_orb, j_orb in product(range(delta.target_shape[0]),range(delta.target_shape[1])): if i_orb != j_orb: @@ -550,7 +543,7 @@ def make_positive_definite(G): name = self.convert_ftps[name.split('_')[0]] solver_eal = solver_eal.real # remove off-diagonal terms - if self.general_params['diag_delta']: + if self.solver_params['diag_delta']: solver_eal = np.diag(np.diag(solver_eal)) h_loc.Fill(name, solver_eal) @@ -561,10 +554,10 @@ def make_positive_definite(G): # so for debugging it is done here again # store solver to h5 archive if self.general_params['store_solver'] and mpi.is_master_node(): - archive = HDFArchive(self.general_params['jobname']+'/'+self.general_params['seedname']+'.h5', 'a') - archive['DMFT_input/solver'].create_group('it_-1') - archive['DMFT_input/solver/it_-1']['Delta'] = self.Delta_freq_solver - archive['DMFT_input/solver/it_-1']['S_'+str(self.icrsh)] = self.triqs_solver + with HDFArchive(self.general_params['jobname']+'/'+self.general_params['seedname']+'.h5', 'a') as archive: + archive['DMFT_input/solver'].create_group('it_-1') + archive['DMFT_input/solver/it_-1']['Delta'] = self.Delta_freq_solver + archive['DMFT_input/solver/it_-1']['S_'+str(self.icrsh)] = self.triqs_solver # Solve the impurity problem for icrsh shell # ************************************* @@ -591,7 +584,7 @@ def make_positive_definite(G): # Solve the impurity problem for icrsh shell # ************************************* - self.triqs_solver.solve(h_int=self.h_int, **{ **self.solver_params, **random_seed }) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # ************************************* # call postprocessing @@ -601,14 +594,13 @@ def make_positive_definite(G): # fill G0_freq from sum_k to solver self.triqs_solver.G0_iw << self.G0_freq - if self.general_params['h_int_type'] == 'dynamic': for b1, b2 in product(self.sum_k.gf_struct_solver_dict[self.icrsh].keys(), repeat=2): self.triqs_solver.D0_iw[b1+"|"+b2] << self.U_iw[self.icrsh] # Solve the impurity problem for icrsh shell # ************************************* - self.triqs_solver.solve(h_int=self.h_int, **{ **self.solver_params, **random_seed }) + self.triqs_solver.solve(h_int=self.h_int, **self.triqs_solver_params) # ************************************* # call postprocessing @@ -626,24 +618,32 @@ def _create_cthyb_solver(self): ''' from triqs_cthyb.solver import Solver as cthyb_solver + # Separately stores all params that go into solve() call of solver + self.triqs_solver_params = {} + keys_to_pass = ('imag_threshold', 'length_cycle', 'max_time', 'measure_density_matrix', + 'measure_G_l', 'measure_pert_order', 'move_double', 'move_shift', + 'n_warmup_cycles', 'loc_n_max', 'loc_n_min', 'off_diag_threshold', 'perform_tail_fit') + for key in keys_to_pass: + self.triqs_solver_params[key] = self.solver_params[key] + # Calculates number of sweeps per rank - self.solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) + self.triqs_solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) # Renames measure chi param - self.solver_params['measure_O_tau_min_ins'] = self.solver_params['measure_chi_insertions'] + self.triqs_solver_params['measure_O_tau_min_ins'] = self.solver_params['measure_chi_insertions'] # use_norm_as_weight also required to measure the density matrix - if self.solver_params['measure_density_matrix']: - self.solver_params['use_norm_as_weight'] = True + self.triqs_solver_params['use_norm_as_weight'] = self.triqs_solver_params['measure_density_matrix'] - # Separately store all params that go into solve() call of solver - self.triqs_solver_params = deepcopy(self.solver_params) - for key in ('type', 'idx_impurities', 'delta_interface', 'fit_max_moment', - 'fit_max_n', 'fit_max_w', 'fit_min_n', 'fit_min_w', 'legendre_fit', - 'measure_chi', 'n_l', 'n_cycles_tot', 'measure_chi_insertions'): - del self.triqs_solver_params[key] - mpi.report(f'Imp {self.icrsh}: passing parameters {self.triqs_solver_params.keys()} ' - 'directly into cthyb_solver.solve()') + if self.triqs_solver_params['perform_tail_fit']: + for key in ('fit_max_moment', 'fit_max_n', 'fit_max_w', 'fit_min_n', 'fit_min_w'): + self.triqs_solver_params[key] = self.solver_params[key] + + # Resolve defaults for n_loc_min and n_loc_max + if self.triqs_solver_params['loc_n_min'] is None: + self.triqs_solver_params['loc_n_min'] = 0 + if self.triqs_solver_params['loc_n_max'] is None: + self.triqs_solver_params['loc_n_max'] = 100000 # large number gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh] # Construct the triqs_solver instances @@ -656,8 +656,6 @@ def _create_cthyb_solver(self): n_iw=self.general_params['n_iw'], n_tau=self.general_params['n_tau'], delta_interface=self.solver_params['delta_interface']) - - return triqs_solver def _create_ctint_solver(self): @@ -666,13 +664,18 @@ def _create_ctint_solver(self): ''' from triqs_ctint import Solver as ctint_solver + # Separately stores all params that go into solve() call of solver + self.triqs_solver_params = {} + keys_to_pass = ('length_cycle', 'max_time', 'measure_pert_order', 'move_double', 'n_warmup_cycles') + for key in keys_to_pass: + self.triqs_solver_params[key] = self.solver_params[key] + # Calculates number of sweeps per rank - self.solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) + self.triqs_solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh] - # TODO: isn't h_int_type a list? - if self.general_params['h_int_type'] == 'dynamic': + if self.general_params['h_int_type'][self.icrsh] == 'dynamic': self.U_iw = None if mpi.is_master_node(): with HDFArchive(self.general_params['jobname']+'/'+self.general_params['seedname']+'.h5', 'r') as archive: @@ -696,6 +699,14 @@ def _create_hubbardI_solver(self): ''' from triqs_hubbardI import Solver as hubbardI_solver + # Separately stores all params that go into solve() call of solver + # All params need to be renamed + self.triqs_solver_params = {} + self.triqs_solver_params['calc_gtau'] = self.solver_params['measure_G_tau'] + self.triqs_solver_params['calc_gw'] = True + self.triqs_solver_params['calc_gl'] = self.solver_params['measure_G_l'] + self.triqs_solver_params['calc_dm'] = self.solver_params['measure_density_matrix'] + gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh] # Construct the triqs_solver instances if self.solver_params['measure_G_l']: @@ -727,8 +738,13 @@ def _create_hartree_solver(self): ''' from triqs_hartree_fock import ImpuritySolver as hartree_solver - gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh] + # Separately stores all params that go into solve() call of solver + self.triqs_solver_params = {} + keys_to_pass = ('method', 'one_shot', 'tol', 'with_fock') + for key in keys_to_pass: + self.triqs_solver_params[key] = self.solver_params[key] + gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh] # Construct the triqs_solver instances # Always initialize the solver with dc_U and dc_J equal to U and J and let the _interface_hartree_dc function # take care of changing the parameters @@ -804,10 +820,29 @@ def _create_ctseg_solver(self): ''' from triqs_ctseg import Solver as ctseg_solver + # Separately stores all params that go into solve() call of solver + self.triqs_solver_params = {} + keys_to_pass = ('improved_estimator', 'length_cycle', 'max_time', 'measure_pert_order', 'n_warmup_cycles') + for key in keys_to_pass: + self.triqs_solver_params[key] = self.solver_params[key] + # Calculates number of sweeps per rank - self.solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) + self.triqs_solver_params['n_cycles'] = int(self.solver_params['n_cycles_tot'] / mpi.size) + + # Rename parameters that are differentin ctseg than cthyb + self.triqs_solver_params['measure_gt'] = self.solver_params['measure_G_tau'] + self.triqs_solver_params['measure_gw'] = self.solver_params['measure_G_iw'] + self.triqs_solver_params['measure_gl'] = self.solver_params['measure_G_l'] + self.triqs_solver_params['measure_hist'] = self.solver_params['measure_pert_order'] + + # Makes sure measure_gw is true if improved estimators are used + if self.triqs_solver_params['improved_estimator']: + self.triqs_solver_params['measure_gt'] = True + self.triqs_solver_params['measure_ft'] = True + else: + self.triqs_solver_params['measure_ft'] = False - if self.general_params['h_int_type'] == 'dynamic': + if self.general_params['h_int_type'][self.icrsh] == 'dynamic': self.U_iw = None if mpi.is_master_node(): with HDFArchive(self.general_params['jobname']+'/'+self.general_params['seedname']+'.h5', 'r') as archive: @@ -836,6 +871,10 @@ def _create_ftps_solver(self): ''' import forktps as ftps + # TODO: add triqs_solver_params for ftps as well to make it analogous to other similars + # Not necessary but nicer. For now, just keep an empty dummy dict + self.triqs_solver_params = {} + # convert self.deg_orbs_ftps to mapping and solver-friendly list if not self.sum_k.corr_shells[self.icrsh]['SO']: # mapping dictionary @@ -921,13 +960,13 @@ def set_Gs_from_G_l(): self.G_time << self.triqs_solver.G_tau self.sum_k.symm_deg_gf(self.G_time, ish=self.icrsh) - if self.general_params['legendre_fit']: + if self.solver_params['legendre_fit']: self.G_time_orig << self.triqs_solver.G_tau # run the filter self.G_l << legendre_filter.apply(self.G_time, self.solver_params['n_l']) # get G_time, G_freq, Sigma_freq from G_l set_Gs_from_G_l() - elif self.solver_params['perform_tail_fit'] and not self.general_params['legendre_fit']: + elif self.solver_params['perform_tail_fit'] and not self.solver_params['legendre_fit']: # if tailfit has been used replace Sigma with the tail fitted Sigma from cthyb self.Sigma_freq << self.triqs_solver.Sigma_iw self.sum_k.symm_deg_gf(self.Sigma_freq, ish=self.icrsh) @@ -977,7 +1016,7 @@ def _ctint_postprocessing(self): self.G_time << Fourier(self.G_freq) # TODO: probably not needed/sensible - # if self.general_params['legendre_fit']: + # if self.solver_params['legendre_fit']: # self.G_freq_orig << self.triqs_solver.G_iw # # run the filter # self.G_l << legendre_filter.apply(self.G_time, self.solver_params['n_l']) @@ -1171,7 +1210,7 @@ def set_Gs_from_G_l(): self.G_l << self.triqs_solver.G_l # get G_time, G_freq, Sigma_freq from G_l set_Gs_from_G_l() - elif self.general_params['legendre_fit']: + elif self.solver_params['legendre_fit']: self.G_time_orig << self.triqs_solver.G_tau self.G_l << legendre_filter.apply(self.G_time, self.solver_params['n_l']) # get G_time, G_freq, Sigma_freq from G_l diff --git a/python/solid_dmft/io_tools/default.toml b/python/solid_dmft/io_tools/default.toml index 2f5d0572..835e5b2d 100644 --- a/python/solid_dmft/io_tools/default.toml +++ b/python/solid_dmft/io_tools/default.toml @@ -11,7 +11,6 @@ csc = false dc = true dc_dmft = "" dc_type = "" -diag_delta = false enforce_off_diag = true eta = "" fixed_mu_value = "" @@ -28,8 +27,6 @@ J = "" jobname = "dmft_dir" load_sigma = false load_sigma_iter = -1 -loc_n_min = "" -loc_n_max = "" magmom = "" magnetic = false mu_gap_gb2_threshold = "" @@ -73,7 +70,9 @@ fit_min_w = "" imag_threshold = 1.0e-14 legendre_fit = false length_cycle = "" -max_time = "" +loc_n_min = "" +loc_n_max = "" +max_time = -1 measure_chi_insertions = 100 measure_chi = "" measure_density_matrix = false @@ -92,7 +91,7 @@ random_seed = "" type = "ctint" idx_impurities = "" length_cycle = "" -max_time = "" +max_time = -1 measure_pert_order = false move_double = true n_cycles_tot = "" @@ -105,7 +104,7 @@ idx_impurities = "" improved_estimator = false legendre_fit = false length_cycle = "" -max_time = "" +max_time = -1 measure_G_iw = false measure_G_l = false measure_G_tau = true @@ -129,6 +128,7 @@ type = "ftps" idx_impurities = "" bath_fit = "" calc_me = true +diag_delta = false dmrg_maxm = 100 dmrg_maxmB = 100 dmrg_maxmI = 100 diff --git a/python/solid_dmft/io_tools/dict_to_h5.py b/python/solid_dmft/io_tools/dict_to_h5.py new file mode 100644 index 00000000..e91307e0 --- /dev/null +++ b/python/solid_dmft/io_tools/dict_to_h5.py @@ -0,0 +1,23 @@ +from copy import deepcopy + +def _iteratively_replace_none(to_write, replace_from, replace_with): + """ Limitation: can only replace None with a string, or a string with None. """ + # First two checks needed because comparison to triqs many-body operator fails + if (isinstance(to_write, str) or to_write is None) and to_write == replace_from: + return replace_with + + if isinstance(to_write, dict): + for key, value in to_write.items(): + to_write[key] = _iteratively_replace_none(value, replace_from, replace_with) + elif isinstance(to_write, list): + for i, value in enumerate(to_write): + to_write[i] = _iteratively_replace_none(value, replace_from, replace_with) + + return to_write + +def prep_params_for_h5(dict_to_write): + return _iteratively_replace_none(deepcopy(dict_to_write), None, 'none') + +# Not sure if the reverse route is actually needed +def prep_params_from_h5(dict_to_read): + return _iteratively_replace_none(deepcopy(dict_to_read), 'none', None) diff --git a/python/solid_dmft/io_tools/verify_input_params.py b/python/solid_dmft/io_tools/verify_input_params.py index e04b9e53..4eb61086 100644 --- a/python/solid_dmft/io_tools/verify_input_params.py +++ b/python/solid_dmft/io_tools/verify_input_params.py @@ -67,32 +67,4 @@ def verify_input_params(params: FullConfig) -> None: def manual_changes_input_params(params: FullConfig) -> None: """ Necessary workarounds for some of the parameters. """ - for entry in params['solver']: - # Calculates the number of solver cycles per rank - # if entry['type'] in ('cthyb', 'ctint', 'ctseg'): - # entry['n_cycles'] = entry['n_cycles_tot'] // mpi.size - # del entry['n_cycles_tot'] - - # Some parameters have different names for ctseg - if entry['type'] == 'ctseg': - entry['measure_gt'] = entry['measure_G_tau'] - del entry['measure_G_tau'] - - entry['measure_gw'] = entry['measure_G_iw'] - del entry['measure_G_iw'] - - # Makes sure measure_gw is true if improved estimators are used - if entry['improved_estimator']: - entry['measure_gt'] = True - entry['measure_ft'] = True - else: - entry['measure_ft'] = False - del entry['improved_estimator'] - - entry['measure_gl'] = entry['measure_G_l'] - del entry['measure_G_l'] - - entry['measure_hist'] = entry['measure_pert_order'] - del entry['measure_pert_order'] - return diff --git a/test/python/CMakeLists.txt b/test/python/CMakeLists.txt index 2c7d6a8a..4ef1020f 100644 --- a/test/python/CMakeLists.txt +++ b/test/python/CMakeLists.txt @@ -1,6 +1,7 @@ # all pytest unittests set (all_pytests test_afm_mapping + test_dict_to_h5 test_interaction_hamiltonian test_manipulate_chemical_potential.py test_observables.py diff --git a/test/python/svo_cthyb_basic_tf/dmft_config.toml b/test/python/svo_cthyb_basic_tf/dmft_config.toml index 2e575d9c..ea95fc26 100644 --- a/test/python/svo_cthyb_basic_tf/dmft_config.toml +++ b/test/python/svo_cthyb_basic_tf/dmft_config.toml @@ -2,7 +2,7 @@ seedname = "inp" jobname = "out" enforce_off_diag = false -block_threshold= 0.001 +block_threshold = 0.001 mu_initial_guess = -0.027041 prec_mu = 0.001 diff --git a/test/python/test_dict_to_h5.py b/test/python/test_dict_to_h5.py new file mode 100644 index 00000000..df408af1 --- /dev/null +++ b/test/python/test_dict_to_h5.py @@ -0,0 +1,13 @@ +from solid_dmft.io_tools import dict_to_h5 + +def test_prep_params_for_h5(): + inp = {'a': None, 'b': {'c': None, 'd': 'e'}, 'f': [None, 'g']} + expected = {'a': 'none', 'b': {'c': 'none', 'd': 'e'}, 'f': ['none', 'g']} + + assert dict_to_h5.prep_params_for_h5(inp) == expected + +def test_prep_params_from_h5(): + inp = {'a': 'none', 'b': {'c': 'none', 'd': 'e'}, 'f': ['none', 'g']} + expected = {'a': None, 'b': {'c': None, 'd': 'e'}, 'f': [None, 'g']} + + assert dict_to_h5.prep_params_from_h5(inp) == expected