update ga code

openms/gwf/__init__.py

@@ -0,0 +1,30 @@
+#
+# @ 2023. Triad National Security, LLC. All rights reserved.
+#
+# This program was produced under U.S. Government contract 89233218CNA000001
+# for Los Alamos National Laboratory (LANL), which is operated by Triad
+# National Security, LLC for the U.S. Department of Energy/National Nuclear
+# Security Administration. All rights in the program are reserved by Triad
+# National Security, LLC, and the U.S. Department of Energy/National Nuclear
+# Security Administration. The Government is granted for itself and others acting
+# on its behalf a nonexclusive, paid-up, irrevocable worldwide license in this
+# material to reproduce, prepare derivative works, distribute copies to the
+# public, perform publicly and display publicly, and to permit others to do so.
+#
+# Author: Yu Zhang <[email protected]>
+#
+
+r"""
+Gutzwiller Waevefunction methods
+================================
+
+This is the collection of GWF solvers for both local, non-local electron-electron correlations
+and electron-boson interactions
+
+"""
+
+# import
+# from openms.gwf import ga_sband
+# from openms.gwf import ga_local
+# from openms.gwf import ga_nonlocal
+# from openms.gwf import ga_dft

openms/gwf/bands.py

@@ -0,0 +1,29 @@
+#
+
+
+def g_fermi(ismear=0, verbose=0):
+    r"""bare band fermi energy
+    TODO:
+    """
+    if ismear == 0:
+        fermi = 0.0
+        # gutz_fermi_tetra_w2k()
+    elif ismear == 1:
+        fermi = 0.0
+        # get_ef_fun()
+    else:
+        # give error
+        raise ValueError("unsupported ismear type")
+    # adjust_efermi_bandgap(verbose)
+
+
+class band_structure(object):
+
+    def __init__(self, *args, **kwargs):
+        self.n_frozen = None  #  number of forzen orbitals
+        self.nelec_frozen = None  #  number of forzen electrons
+        self.ne = None #  
+        self.ek = None # band eigen values (nband, nk, nspin)
+        self.coeff = None # spin-up/down weights
+        self.efermi = 0.0 # fermi energy
+        self.hk0 = None

openms/gwf/ga_dft.py

@@ -0,0 +1,101 @@
+from pyscf.pbc.scf import khf
+from pyscf.pbc.dft import rks
+from pyscf.pbc.dft import krks
+
+
+# suggestion from chatgpt (not usefull at all)
+r"""
+
+This pseudocode outlines a basic structure for implementing LDA+Gutzwiller in PySCF. Here's a breakdown of the key components:
+
+Gutzwiller-related functions: These include the projector, energy calculation, potential calculation, parameter optimization,
+and renormalization factors. These functions need to be implemented in detail.
+LDAGutzwiller class: This subclasses PySCF's RKS (Restricted Kohn-Sham) class to include Gutzwiller corrections. The key modifications are:
+
+Initialization with Hubbard U and Hund's J parameters.
+Modified get_veff method to include the Gutzwiller potential.
+Modified energy_elec method to include the Gutzwiller energy contribution.
+Modified scf method to optimize Gutzwiller parameters in each SCF cycle.
+"""
+
+
+# Define Gutzwiller-related functions
+def gutzwiller_projector(wavefunction, g_params):
+    # Apply Gutzwiller projection to wavefunctions
+    # ...
+    pass
+
+
+def gutzwiller_energy(dm, g_params, U, J):
+    # Calculate Gutzwiller energy contribution
+    # ...
+    pass
+
+
+def gutzwiller_potential(dm, g_params):
+    # Calculate Gutzwiller potential
+    # ...
+    pass
+
+
+def optimize_gutzwiller_params(wavefunction, dm, U, J):
+    # Optimize Gutzwiller parameters
+    # ...
+    pass
+
+
+def renormalization_factors(g_params, occupations):
+    # Calculate renormalization factors
+    # ...
+    pass
+
+
+def gutzwiller_corrected_hamiltonian(h_core, g_potential):
+    # Construct Gutzwiller-corrected Hamiltonian
+    # ...
+    pass
+
+
+class GRKS(krks.KRKS):
+    r"""
+    LDA+GUTZWILLER
+    """
+
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+
+    def get_veff(self, mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1):
+        # Get standard LDA effective potential
+        veff = super().get_veff(mol, dm, dm_last, vhf_last, hermi)
+
+        # Add Gutzwiller potential
+        v_gutzwiller = gutzwiller_potential(dm, self.g_params)
+        veff += v_gutzwiller
+
+        return veff
+
+    def energy_elec(self, dm=None, h1e=None, vhf=None):
+        # Get standard LDA electronic energy
+        e_tot, e_coul = super().energy_elec(dm, h1e, vhf)
+
+        # Add Gutzwiller energy contribution
+        e_gutzwiller = gutzwiller_energy(dm, self.g_params, self.U, self.J)
+        e_tot += e_gutzwiller
+
+        return e_tot, e_coul
+
+    def scf(self, *args, **kwargs):
+        for cycle in range(self.max_cycle):
+            # Perform standard SCF iteration
+            dm = super().scf(*args, **kwargs)
+
+            # Optimize Gutzwiller parameters
+            self.g_params = optimize_gutzwiller_params(
+                self.mo_coeff, dm, self.U, self.J
+            )
+
+            # Check for convergence
+            if self.converged:
+                break
+
+        return dm