From b9e141c7c882ef13e2aab95bd76f0034f3830ceb Mon Sep 17 00:00:00 2001 From: julpe Date: Fri, 27 Mar 2026 09:43:09 +0100 Subject: [PATCH] Fixed the Pulay mixing, it should be more stable now. --- moldga/__init__.py | 6 + moldga/brillouin_zone.py | 9 +- moldga/bubble_gen.py | 6 + moldga/config.py | 6 + moldga/config_parser.py | 10 +- moldga/dga_io.py | 6 + moldga/dga_logger.py | 6 + moldga/dga_main.py | 6 + moldga/eliashberg_solver.py | 6 + moldga/four_point.py | 6 + moldga/gap_function.py | 6 + moldga/greens_function.py | 6 + moldga/hamiltonian.py | 6 + moldga/interaction.py | 6 + moldga/lambda_correction.py | 6 + moldga/local_four_point.py | 6 + moldga/local_n_point.py | 6 + moldga/local_sde.py | 6 + moldga/matsubara_frequencies.py | 6 + moldga/mpi_distributor.py | 6 + moldga/n_point_base.py | 61 ++++-- moldga/nonlocal_sde.py | 31 ++- moldga/plotting.py | 6 + moldga/self_energy.py | 6 + moldga/symmetrize_new.py | 6 + moldga/w2dyn_aux.py | 6 + tests/__init__.py | 6 + tests/conftest.py | 6 + tests/test_brillouin_zone.py | 6 + tests/test_eliashberg_end_to_end.py | 6 + tests/test_four_point.py | 6 + tests/test_greens_function.py | 6 + tests/test_hamiltonian.py | 6 + tests/test_interaction.py | 6 + tests/test_lambda_correction.py | 6 + tests/test_local_four_point.py | 6 + tests/test_local_n_point.py | 6 + tests/test_local_sde_end_to_end.py | 6 + tests/test_n_point_base.py | 6 + tests/test_nonlocal_sde_end_to_end.py | 21 +- tests/test_pulay_mixing.py | 300 ++++++++++++++++++++++++++ tests/test_self_energy.py | 6 + tests/test_symmetrize.py | 6 + 43 files changed, 617 insertions(+), 37 deletions(-) create mode 100644 tests/test_pulay_mixing.py diff --git a/moldga/__init__.py b/moldga/__init__.py index e69de29..f4d5a35 100644 --- a/moldga/__init__.py +++ b/moldga/__init__.py @@ -0,0 +1,6 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + diff --git a/moldga/brillouin_zone.py b/moldga/brillouin_zone.py index b6e1f9f..20876e2 100644 --- a/moldga/brillouin_zone.py +++ b/moldga/brillouin_zone.py @@ -1,6 +1,11 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + """ -Module to handle operations within the (irreducible) Brillouin zone. Copied over from Paul Worm's code. -Only modified the constant arrays and made enums out of them for type hinting. +Module to handle operations within the (irreducible) Brillouin zone. Heavily inspired by DGApy. """ import warnings diff --git a/moldga/bubble_gen.py b/moldga/bubble_gen.py index cc7793a..4ced85b 100644 --- a/moldga/bubble_gen.py +++ b/moldga/bubble_gen.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import moldga.config as config from moldga.greens_function import GreensFunction from moldga.matsubara_frequencies import MFHelper diff --git a/moldga/config.py b/moldga/config.py index 097ed28..e90e808 100644 --- a/moldga/config.py +++ b/moldga/config.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import numpy as np import moldga.brillouin_zone as bz diff --git a/moldga/config_parser.py b/moldga/config_parser.py index 303b61e..8f17432 100644 --- a/moldga/config_parser.py +++ b/moldga/config_parser.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import argparse import os @@ -207,8 +213,8 @@ def _build_self_energy_interpolation_config(self, conf_file): section = conf_file["self_energy_interpolation"] conf.do_interpolation = self._try_parse(section, "do_interpolation", False) - conf.beta_target = self._try_parse(section, "beta_target", 1.0) - conf.niv_target = self._try_parse(section, "niv_target", 10) + conf.beta_target = self._try_parse(section, "target_beta", 1.0) + conf.niv_target = self._try_parse(section, "target_niv", 10) return conf diff --git a/moldga/dga_io.py b/moldga/dga_io.py index 1a95814..6863f9e 100644 --- a/moldga/dga_io.py +++ b/moldga/dga_io.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os import moldga.brillouin_zone as bz diff --git a/moldga/dga_logger.py b/moldga/dga_logger.py index bf1cf41..bdd1b92 100644 --- a/moldga/dga_logger.py +++ b/moldga/dga_logger.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import logging import os from datetime import datetime diff --git a/moldga/dga_main.py b/moldga/dga_main.py index 725e5a7..f2e6076 100644 --- a/moldga/dga_main.py +++ b/moldga/dga_main.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import itertools as it import logging import os diff --git a/moldga/eliashberg_solver.py b/moldga/eliashberg_solver.py index f8bff3d..600da42 100644 --- a/moldga/eliashberg_solver.py +++ b/moldga/eliashberg_solver.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os from typing import Tuple diff --git a/moldga/four_point.py b/moldga/four_point.py index c54b7f2..41de18f 100644 --- a/moldga/four_point.py +++ b/moldga/four_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import gc from copy import deepcopy diff --git a/moldga/gap_function.py b/moldga/gap_function.py index 956162f..a64e343 100644 --- a/moldga/gap_function.py +++ b/moldga/gap_function.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import numpy as np from moldga.brillouin_zone import KGrid diff --git a/moldga/greens_function.py b/moldga/greens_function.py index 5fb32bb..db80498 100644 --- a/moldga/greens_function.py +++ b/moldga/greens_function.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from copy import deepcopy import numpy as np diff --git a/moldga/hamiltonian.py b/moldga/hamiltonian.py index 69b202e..7505107 100644 --- a/moldga/hamiltonian.py +++ b/moldga/hamiltonian.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import itertools as it import logging diff --git a/moldga/interaction.py b/moldga/interaction.py index aab3c58..7bfd372 100644 --- a/moldga/interaction.py +++ b/moldga/interaction.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from copy import deepcopy import numpy as np diff --git a/moldga/lambda_correction.py b/moldga/lambda_correction.py index 75f2569..2f22c3f 100644 --- a/moldga/lambda_correction.py +++ b/moldga/lambda_correction.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import numpy as np from moldga import config diff --git a/moldga/local_four_point.py b/moldga/local_four_point.py index bc963c2..dd8d8fa 100644 --- a/moldga/local_four_point.py +++ b/moldga/local_four_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from moldga.interaction import LocalInteraction, Interaction from moldga.local_n_point import LocalNPoint from moldga.matsubara_frequencies import MFHelper diff --git a/moldga/local_n_point.py b/moldga/local_n_point.py index b92073b..32d45ad 100644 --- a/moldga/local_n_point.py +++ b/moldga/local_n_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import itertools import os diff --git a/moldga/local_sde.py b/moldga/local_sde.py index 340ef8a..40be5f4 100644 --- a/moldga/local_sde.py +++ b/moldga/local_sde.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import moldga.config as config from moldga.bubble_gen import BubbleGenerator from moldga.greens_function import GreensFunction diff --git a/moldga/matsubara_frequencies.py b/moldga/matsubara_frequencies.py index 909dc2b..56c6637 100644 --- a/moldga/matsubara_frequencies.py +++ b/moldga/matsubara_frequencies.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from enum import Enum import numpy as np diff --git a/moldga/mpi_distributor.py b/moldga/mpi_distributor.py index 5a570c2..5575966 100644 --- a/moldga/mpi_distributor.py +++ b/moldga/mpi_distributor.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import gc import os diff --git a/moldga/n_point_base.py b/moldga/n_point_base.py index 17e268e..77fabd9 100644 --- a/moldga/n_point_base.py +++ b/moldga/n_point_base.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import gc from abc import ABC from copy import deepcopy @@ -490,29 +496,54 @@ def _map_to_full_bz(self, k_grid: KGrid, num_orbital_dimensions: int, nq: tuple if np.allclose(u_ref, identity): continue - uc_ref = u_ref.conj() idx = np.array(indices) - if num_orbital_dimensions == 2: - if path_2 is None: - path_2 = np.einsum_path( - "ap,bq,kpq...->kab...", u_ref, uc_ref, self.mat[idx], optimize="optimal" - )[0] - self.mat[idx] = np.einsum("ap,bq,kpq...->kab...", u_ref, uc_ref, self.mat[idx], optimize=path_2) - elif num_orbital_dimensions == 4: - if path_4 is None: - path_4 = np.einsum_path( + def _is_permutation_matrix(u: np.ndarray) -> bool: + return ( + np.allclose(np.abs(u), np.abs(u).astype(int)) # only 0s and 1s + and np.allclose(u.sum(axis=0), 1) # one 1 per column + and np.allclose(u.sum(axis=1), 1) # one 1 per row + ) + + if _is_permutation_matrix(u_ref): + # For real permutation matrices, U @ M @ U^T is just index reordering. + # Find the permutation: perm[i] = j means orbital i gets content from orbital j. + perm = np.argmax(u_ref.real, axis=1) + + if num_orbital_dimensions == 2: + self.mat[idx] = self.mat[idx][:, perm, ...][:, :, perm, ...] + elif num_orbital_dimensions == 4: + self.mat[idx] = self.mat[idx][:, perm, ...][:, :, perm, ...][:, :, :, perm, ...][ + :, :, :, :, perm, ... + ] + else: + uc_ref = u_ref.conj() + if num_orbital_dimensions == 2: + if path_2 is None: + path_2 = np.einsum_path( + "ap,bq,kpq...->kab...", u_ref, uc_ref, self.mat[idx], optimize="optimal" + )[0] + self.mat[idx] = np.einsum("ap,bq,kpq...->kab...", u_ref, uc_ref, self.mat[idx], optimize=path_2) + elif num_orbital_dimensions == 4: + if path_4 is None: + path_4 = np.einsum_path( + "ap,bq,cr,ds,kpqrs...->kabcd...", + u_ref, + uc_ref, + u_ref, + uc_ref, + self.mat[idx], + optimize="optimal", + )[0] + self.mat[idx] = np.einsum( "ap,bq,cr,ds,kpqrs...->kabcd...", u_ref, uc_ref, u_ref, uc_ref, self.mat[idx], - optimize="optimal", - )[0] - self.mat[idx] = np.einsum( - "ap,bq,cr,ds,kpqrs...->kabcd...", u_ref, uc_ref, u_ref, uc_ref, self.mat[idx], optimize=path_4 - ) + optimize=path_4, + ) self.mat = self.mat.reshape((np.prod(self.nq), *self.original_shape[1:])) self.update_original_shape() diff --git a/moldga/nonlocal_sde.py b/moldga/nonlocal_sde.py index 5d36473..83c269c 100644 --- a/moldga/nonlocal_sde.py +++ b/moldga/nonlocal_sde.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import gc import glob import os @@ -704,11 +710,11 @@ def apply_mixing_strategy( is either 'linear' or 'pulay'. """ logger = config.logger - strategy = config.self_consistency.mixing_strategy n_hist = config.self_consistency.mixing_history_length + alpha = config.self_consistency.mixing if ( - strategy == "pulay" + config.self_consistency.mixing_strategy == "pulay" and current_iter > n_hist and config.self_consistency.save_iter and config.output.save_quantities @@ -720,10 +726,10 @@ def apply_mixing_strategy( last_proposals = [sigma_dmft_stacked] + last_results last_results = last_results + [sigma_new.mat] - niv_dmft = sigma_new.current_shape[-1] // 2 + niv = sigma_new.current_shape[-1] // 2 niv_core = config.box.niv_core - last_proposals = [sigma[..., niv_dmft - niv_core : niv_dmft + niv_core] for sigma in last_proposals] - last_results = [sigma[..., niv_dmft - niv_core : niv_dmft + niv_core] for sigma in last_results] + last_proposals = [sigma[..., niv - niv_core : niv + niv_core] for sigma in last_proposals] + last_results = [sigma[..., niv - niv_core : niv + niv_core] for sigma in last_results] logger.info(f"Loaded last {n_hist} self-energies from files.") shape = last_results[-1].shape @@ -749,17 +755,20 @@ def get_result(idx: int): f_matrix[:, i] -= r_matrix[:, i] + # Residual: F(x_n) - x_n, where x_n = last_proposals[-1] = sigma_old (core window) iter_diff = get_result(-1) - get_proposal(-1) f_i[:n_total] = iter_diff.real f_i[n_total:] = iter_diff.imag - update = config.self_consistency.mixing * f_i - fact1 = (r_matrix + config.self_consistency.mixing * f_matrix) @ np.linalg.inv(f_matrix.T @ f_matrix) - update -= fact1 @ (f_matrix.T @ f_i) + # Solve min||F @ c - f_i|| via least squares (more stable than explicit inverse) + coeffs, _, _, _ = np.linalg.lstsq(f_matrix, f_i, rcond=None) + + # Pulay update: x_{n+1} = x_n + alpha*f_i - (R + alpha*F) @ c + update = alpha * f_i - (r_matrix + alpha * f_matrix) @ coeffs update = update[:n_total] + 1j * update[n_total:] - sigma_new.mat[..., niv_dmft - niv_core : niv_dmft + niv_core] = sigma_old.compress_q_dimension().mat[ - ..., niv_dmft - niv_core : niv_dmft + niv_core - ] + update.reshape(shape) + + # Update the new self energy + sigma_new.mat[..., niv - niv_core : niv + niv_core] = get_proposal(-1).reshape(shape) + update.reshape(shape) logger.info( f"Pulay mixing applied with {n_hist} previous iterations and a mixing parameter of {config.self_consistency.mixing}." diff --git a/moldga/plotting.py b/moldga/plotting.py index 99ff156..b2e0bc3 100644 --- a/moldga/plotting.py +++ b/moldga/plotting.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os import numpy as np diff --git a/moldga/self_energy.py b/moldga/self_energy.py index e4b89ae..ac2f098 100644 --- a/moldga/self_energy.py +++ b/moldga/self_energy.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import itertools import itertools as it from copy import deepcopy diff --git a/moldga/symmetrize_new.py b/moldga/symmetrize_new.py index e5440b9..42d37fb 100644 --- a/moldga/symmetrize_new.py +++ b/moldga/symmetrize_new.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import gc import glob import itertools as it diff --git a/moldga/w2dyn_aux.py b/moldga/w2dyn_aux.py index 58c2a7f..ac59aff 100644 --- a/moldga/w2dyn_aux.py +++ b/moldga/w2dyn_aux.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import h5py import moldga.symmetrize_new as sym diff --git a/tests/__init__.py b/tests/__init__.py index e69de29..f4d5a35 100644 --- a/tests/__init__.py +++ b/tests/__init__.py @@ -0,0 +1,6 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + diff --git a/tests/conftest.py b/tests/conftest.py index 2d50e46..b18d127 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import logging import os from unittest.mock import MagicMock diff --git a/tests/test_brillouin_zone.py b/tests/test_brillouin_zone.py index 30364fc..f5ecf94 100644 --- a/tests/test_brillouin_zone.py +++ b/tests/test_brillouin_zone.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import warnings import numpy as np diff --git a/tests/test_eliashberg_end_to_end.py b/tests/test_eliashberg_end_to_end.py index 36a5e48..a45f5c1 100644 --- a/tests/test_eliashberg_end_to_end.py +++ b/tests/test_eliashberg_end_to_end.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os from unittest.mock import patch diff --git a/tests/test_four_point.py b/tests/test_four_point.py index 8340869..2861bb9 100644 --- a/tests/test_four_point.py +++ b/tests/test_four_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from unittest.mock import MagicMock import numpy as np diff --git a/tests/test_greens_function.py b/tests/test_greens_function.py index 4f6e8d5..a5eccb7 100644 --- a/tests/test_greens_function.py +++ b/tests/test_greens_function.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from unittest.mock import MagicMock import numpy as np diff --git a/tests/test_hamiltonian.py b/tests/test_hamiltonian.py index f00ff55..70b0a78 100644 --- a/tests/test_hamiltonian.py +++ b/tests/test_hamiltonian.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os import numpy as np diff --git a/tests/test_interaction.py b/tests/test_interaction.py index c2339d5..69d5f55 100644 --- a/tests/test_interaction.py +++ b/tests/test_interaction.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import pytest import numpy as np from moldga.interaction import Interaction, LocalInteraction, SpinChannel diff --git a/tests/test_lambda_correction.py b/tests/test_lambda_correction.py index 9f9fb8a..0184ef5 100644 --- a/tests/test_lambda_correction.py +++ b/tests/test_lambda_correction.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os from unittest.mock import patch diff --git a/tests/test_local_four_point.py b/tests/test_local_four_point.py index 803e90b..4aa96b9 100644 --- a/tests/test_local_four_point.py +++ b/tests/test_local_four_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from unittest.mock import patch, MagicMock import pytest diff --git a/tests/test_local_n_point.py b/tests/test_local_n_point.py index 3ac742a..5659a26 100644 --- a/tests/test_local_n_point.py +++ b/tests/test_local_n_point.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import itertools from unittest.mock import patch diff --git a/tests/test_local_sde_end_to_end.py b/tests/test_local_sde_end_to_end.py index 336b773..1fa5538 100644 --- a/tests/test_local_sde_end_to_end.py +++ b/tests/test_local_sde_end_to_end.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import logging import os from unittest.mock import patch, MagicMock diff --git a/tests/test_n_point_base.py b/tests/test_n_point_base.py index 79aac13..25e4a89 100644 --- a/tests/test_n_point_base.py +++ b/tests/test_n_point_base.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import os import sys import types diff --git a/tests/test_nonlocal_sde_end_to_end.py b/tests/test_nonlocal_sde_end_to_end.py index eafad23..937ed76 100644 --- a/tests/test_nonlocal_sde_end_to_end.py +++ b/tests/test_nonlocal_sde_end_to_end.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import contextlib import os import types @@ -36,13 +42,12 @@ def create_srvo3_cubic_config(c, f: str): c.box.niv_shell = 0 c.output.save_quantities = False c.output.do_plotting = False - c.lattice.nk = (8, 8, 8) + c.lattice.nk = (12, 12, 12) c.lattice.nq = config.lattice.nk - c.lattice.symmetries = "three_dimensional_cubic" c.lattice.k_grid = bz.KGrid(c.lattice.nk, bz.three_dimensional_cubic_symmetries()) - c.lattice.q_grid = config.lattice.k_grid - c.lattice.k_grid.specify_orbital_basis(3, "t2g") - c.lattice.q_grid.specify_orbital_basis(3, "t2g") + c.lattice.q_grid = c.lattice.k_grid + c.lattice.symmetries = "three_dimensional_cubic" + c.lattice.orbital_basis = "t2g" c.lattice.type = "from_wannier90" c.lattice.interaction_type = "kanamori_from_dmft" c.lattice.er_input = f"{f}/wan_hr.dat" @@ -188,13 +193,15 @@ def test_calculates_srvo3_correctly(setup_srvo3_cubic): sigma_dga_cubic = nonlocal_sde.calculate_self_energy_q(comm_mock, u_loc, v_nonloc, s_dmft, s_loc) niv = sigma_dga_cubic.current_shape[-1] // 2 - s_cubic = sigma_dga_cubic.compress_q_dimension().mat.reshape(8, 8, 8, 3, 3, 2 * niv) # (nkx, nky, nkz, nb, nb, niv) + s_cubic = sigma_dga_cubic.compress_q_dimension().mat.reshape( + 12, 12, 12, 3, 3, 2 * niv + ) # (nkx, nky, nkz, nb, nb, niv) s_xy_cub = np.swapaxes(s_cubic, 0, 1) s_xz_cub = np.swapaxes(s_cubic, 0, 2) s_yz_cub = np.swapaxes(s_cubic, 1, 2) - atol = 1e-2 # 8x8x8 grid is too coarse for tight mirror symmetry checks + atol = 1e-6 # X_Y_SYM (kx<->ky): dxy(0) invariant, dxz(1)<->dyz(2) swap assert np.allclose(s_cubic[..., 0, 0, :], s_xy_cub[..., 0, 0, :], atol=atol), "X_Y_SYM dxy failed" diff --git a/tests/test_pulay_mixing.py b/tests/test_pulay_mixing.py new file mode 100644 index 0000000..2369baf --- /dev/null +++ b/tests/test_pulay_mixing.py @@ -0,0 +1,300 @@ +import numpy as np +import pytest +from unittest.mock import patch + +import moldga.config as real_config +from moldga.self_energy import SelfEnergy +from moldga.nonlocal_sde import apply_mixing_strategy + + +# --------------------------------------------------------------------------- +# Module-level setup: patch config.sys.beta so SelfEnergy.__init__ works +# --------------------------------------------------------------------------- + +BETA = 10.0 +NB = 1 +NK = (1, 1, 1) +NIV = 8 +NIV_CORE = 4 + + +@pytest.fixture(autouse=True) +def set_beta(): + """Patches config.sys.beta globally for all tests so SelfEnergy can be constructed.""" + original = real_config.sys.beta + real_config.sys.beta = BETA + yield + real_config.sys.beta = original + + +# --------------------------------------------------------------------------- +# Helpers +# --------------------------------------------------------------------------- + + +def make_sigma(value: complex, nk: tuple = NK, nb: int = NB, niv: int = NIV) -> SelfEnergy: + """Creates a SelfEnergy with constant complex fill value.""" + mat = np.full((*nk, nb, nb, 2 * niv), value, dtype=np.complex64) + return SelfEnergy(mat, nk, full_niv_range=True, has_compressed_q_dimension=True) + + +def make_sigma_mat(value: complex, nk: tuple = NK, nb: int = NB, niv: int = NIV) -> np.ndarray: + """Returns a raw numpy array with the given fill value in the expected shape.""" + nk_tot = int(np.prod(nk)) + return np.full((nk_tot, nb, nb, 2 * niv), value, dtype=np.complex64) + + +def make_config_mock( + strategy: str = "linear", + mixing: float = 0.5, + n_hist: int = 3, + save_iter: bool = True, + save_quantities: bool = True, + niv_core: int = NIV_CORE, + nk_tot: int = 1, + output_path: str = "./", + previous_sc_path: str = "./", +): + """Builds a mock config object for patching moldga.nonlocal_sde.config.""" + from unittest.mock import MagicMock + + cfg = MagicMock() + cfg.self_consistency.mixing_strategy = strategy + cfg.self_consistency.mixing = mixing + cfg.self_consistency.mixing_history_length = n_hist + cfg.self_consistency.save_iter = save_iter + cfg.self_consistency.previous_sc_path = previous_sc_path + cfg.output.save_quantities = save_quantities + cfg.output.output_path = output_path + cfg.box.niv_core = niv_core + cfg.lattice.k_grid.nk_tot = nk_tot + cfg.logger = MagicMock() + return cfg + + +def patch_config(**kwargs): + return patch("moldga.nonlocal_sde.config", make_config_mock(**kwargs)) + + +def run_pulay( + sigma_new: SelfEnergy, + sigma_old: SelfEnergy, + sigma_dmft: SelfEnergy, + file_sigmas: list, + mixing: float = 0.5, + n_hist: int = 3, + niv_core: int = NIV_CORE, + current_iter: int = 10, +) -> SelfEnergy: + nk_tot = int(np.prod(NK)) + with ( + patch_config(strategy="pulay", mixing=mixing, n_hist=n_hist, niv_core=niv_core, nk_tot=nk_tot), + patch("moldga.nonlocal_sde.read_last_n_sigmas_from_files", return_value=file_sigmas), + ): + return apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=current_iter) + + +# --------------------------------------------------------------------------- +# Linear mixing tests +# --------------------------------------------------------------------------- + + +def test_linear_mixing_basic(): + """x_mixed = alpha * x_new + (1 - alpha) * x_old""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(0.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="linear", mixing=0.5): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=1) + + np.testing.assert_allclose(result.mat, 1.0, atol=1e-5) + + +def test_linear_mixing_alpha_zero(): + """alpha=0 should return sigma_old unchanged.""" + sigma_new = make_sigma(5.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="linear", mixing=0.0): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=1) + + np.testing.assert_allclose(result.mat, 1.0, atol=1e-5) + + +def test_linear_mixing_alpha_one(): + """alpha=1 should return sigma_new unchanged.""" + sigma_new = make_sigma(5.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="linear", mixing=1.0): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=1) + + np.testing.assert_allclose(result.mat, 5.0, atol=1e-5) + + +def test_linear_mixing_complex(): + """Linear mixing should work correctly for complex-valued self-energies.""" + sigma_new = make_sigma(2.0 + 2.0j) + sigma_old = make_sigma(0.0 + 0.0j) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="linear", mixing=0.5): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=1) + + np.testing.assert_allclose(result.mat, 1.0 + 1.0j, atol=1e-5) + + +def test_linear_mixing_returns_self_energy_instance(): + """Linear mixing must return a SelfEnergy instance.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="linear", mixing=0.5): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=1) + + assert isinstance(result, SelfEnergy) + + +# --------------------------------------------------------------------------- +# Pulay fallback tests +# --------------------------------------------------------------------------- + + +def test_pulay_falls_back_to_linear_when_iter_too_small(): + """Pulay mixing must fall back to linear if current_iter <= n_hist.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(0.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="pulay", mixing=0.5, n_hist=5): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=3) + + np.testing.assert_allclose(result.mat, 1.0, atol=1e-5) + + +def test_pulay_falls_back_to_linear_when_save_iter_false(): + """Pulay requires save_iter=True, otherwise falls back to linear.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(0.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="pulay", mixing=0.5, n_hist=3, save_iter=False): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=10) + + np.testing.assert_allclose(result.mat, 1.0, atol=1e-5) + + +def test_pulay_falls_back_to_linear_when_save_quantities_false(): + """Pulay requires save_quantities=True, otherwise falls back to linear.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(0.0) + sigma_dmft = make_sigma(0.0) + + with patch_config(strategy="pulay", mixing=0.5, n_hist=3, save_quantities=False): + result = apply_mixing_strategy(sigma_new, sigma_old, sigma_dmft, current_iter=10) + + np.testing.assert_allclose(result.mat, 1.0, atol=1e-5) + + +# --------------------------------------------------------------------------- +# Pulay mixing tests +# --------------------------------------------------------------------------- + + +def test_pulay_returns_self_energy_instance(): + """Pulay mixing must return a SelfEnergy instance.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(float(v)) for v in [0.5, 0.8, 1.0]] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + assert isinstance(result, SelfEnergy) + + +def test_pulay_converged_fixed_point(): + """If all sigmas are identical, Pulay mixing must return the same sigma in the core window.""" + value = 3.0 + 1.0j + sigma_new = make_sigma(value) + sigma_old = make_sigma(value) + sigma_dmft = make_sigma(value) + file_sigmas = [make_sigma_mat(value) for _ in range(3)] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + niv_dmft = sigma_new.mat.shape[-1] // 2 + np.testing.assert_allclose( + result.mat[..., niv_dmft - NIV_CORE : niv_dmft + NIV_CORE], + np.full_like(result.mat[..., niv_dmft - NIV_CORE : niv_dmft + NIV_CORE], value), + atol=1e-4, + ) + + +def test_pulay_returns_same_object_as_sigma_new(): + """Pulay mixing writes into sigma_new directly and returns it.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(float(v)) for v in [0.5, 0.8, 1.0]] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + assert result is sigma_new + + +def test_pulay_does_not_mutate_sigma_old(): + """apply_mixing_strategy must not corrupt sigma_old.mat.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(1.0) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(float(v)) for v in [0.5, 0.8, 1.0]] + + original_mat = sigma_old.mat.copy() + run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + np.testing.assert_array_equal(sigma_old.mat, original_mat) + + +def test_pulay_tails_come_from_sigma_new(): + """Frequencies outside the core window must be taken from sigma_new, not sigma_old.""" + sigma_new = make_sigma(2.0) + sigma_old = make_sigma(99.0) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(2.0) for _ in range(3)] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + niv_dmft = sigma_new.mat.shape[-1] // 2 + np.testing.assert_allclose(result.mat[..., : niv_dmft - NIV_CORE], 2.0, atol=1e-5) + np.testing.assert_allclose(result.mat[..., niv_dmft + NIV_CORE :], 2.0, atol=1e-5) + + +def test_pulay_result_shape_matches_sigma_new(): + """The result must have the same shape as sigma_new.mat.""" + sigma_new = make_sigma(1.0) + sigma_old = make_sigma(0.5) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(float(v)) for v in [0.3, 0.4, 0.5]] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + assert result.mat.shape == sigma_new.mat.shape + + +def test_pulay_core_is_finite(): + """The core window of the Pulay result must contain only finite values.""" + sigma_new = make_sigma(1.5 + 0.5j) + sigma_old = make_sigma(1.0 + 0.3j) + sigma_dmft = make_sigma(0.0) + file_sigmas = [make_sigma_mat(complex(0.5 + 0.1j * i)) for i in range(3)] + + result = run_pulay(sigma_new, sigma_old, sigma_dmft, file_sigmas) + + niv_dmft = result.mat.shape[-1] // 2 + core = result.mat[..., niv_dmft - NIV_CORE : niv_dmft + NIV_CORE] + assert np.all(np.isfinite(core)) diff --git a/tests/test_self_energy.py b/tests/test_self_energy.py index 87cb610..b4876af 100644 --- a/tests/test_self_energy.py +++ b/tests/test_self_energy.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + from unittest.mock import MagicMock import numpy as np diff --git a/tests/test_symmetrize.py b/tests/test_symmetrize.py index ca3b985..0bd3fc5 100644 --- a/tests/test_symmetrize.py +++ b/tests/test_symmetrize.py @@ -1,3 +1,9 @@ +# SPDX-FileCopyrightText: 2025-2026 Julian Peil +# SPDX-License-Identifier: MIT +# +# moLDGA — Multi-Orbital Ladder Dynamical Vertex Approximation (LDGA) & +# Eliashberg Equation Solver for Strongly Correlated Electron Systems + import pytest from moldga.symmetrize_new import *