pyseqm interface fixes

examples/pyseqm_interface/single_node/input.in

@@ -5,12 +5,12 @@ Threshold= 1.0E-5                                              # Currently ignor
 CoordsFile= nanostar.pdb                                       # Geometry filename.
 GraphThreshold= 0.0005                                         # Density matrix graph threshold (frobenius norm of off-diagonal blocks).
 GraphThresholdInit= 0.04                                       # For initial overlap matrix graph. Ignored if InitGraphType= Geom.
-MaxDeg= 500                                                    # Max graph degree. Must be large enough to fit the largest core+halo.
+MaxDeg= 900                                                    # Max graph degree. Must be large enough to fit the largest core+halo.
 Rcut= 3.0                                                      # Radius cutoff for initial graph. Ignored if InitGraphType= OverlapM.
 Alpha= 0.7                                                     # Density matrix linear mixing coeff. Alpha*new + (1-Alpha)*old.
 InitGraphType= Geom                                            # How to generate initial graph. Geom: through distances. OverlapM: through overlap matrix.
 PartitionType= SpectralClustering                              # SpectralClustering MinCut BalancedKmeans.
-NumParts= 4                                                    # Number of cores.
+NumParts= 4                                                  # Number of cores.
 SCFTol= 1.0E-5                                                 # Currently ignored
 Orbitals= {"H":1,"Li":4,"C":4,"N":4,"O":4,"F":4,"Cl":4,"Na":4} # This has to be known for each engine.
 Valency= {"H":1,"Li":1,"C":4,"N":5,"O":6,"F":7,"Cl":7,"Na":1}  # This has to be known for each engine.
@@ -24,5 +24,5 @@ restartSave= False                                             # Save files for
 restartLoad= False                                             # Load files for scf restart.
 numGPU= -1                                                     # If on CPU, no efect. If on GPU and nodes are homogeneous, set to -1. If on GPU and nodes are inhomogeneous (diferent number of GPUs per node), set it to the minimum number of GPUs per node. 
 writeGeom= False                                               # Save core and core+halo geometries.
-Engine= {"Name":"ProxyAPython","InterfaceType":"PySEQM","Path":".../sedacs/proxies/python/","Executable":""} #The name of the external code
+Engine= {"Name":"XXX","InterfaceType":"PySEQM","Path":"/home/maxim/Projects/SEDACS_github/sedacs/proxies/python/","Executable":""} #The name of the external code
 

examples/pyseqm_interface/single_node/main.py

@@ -7,9 +7,17 @@
 import os
 
 # ADD PROXYA_PATH to PYTHONPATH
-sedacs_base = importlib.util.find_spec("sedacs").submodule_search_locations[0]
-target_path = os.path.abspath(os.path.join(sedacs_base, "../../"))
-sys.path.append(target_path)
+# sedacs_base = importlib.util.find_spec("sedacs").submodule_search_locations[0]
+# target_path = os.path.abspath(os.path.join(sedacs_base, "../../"))
+# sys.path.append(target_path)
+
+# ADD PROXYA_PATH to PYTHONPATH
+proxya_path = "/home/maxim/Projects/SEDACS_github/sedacs"
+sys.path.append(proxya_path)
+
+### ADD PATH TO PYSEQM ###
+pyseqm_path = "/home/maxim/Projects/git2/PYSEQM_dev/"
+sys.path.insert(1, pyseqm_path)
 
 DTYPE = torch.float64
 torch.set_default_dtype(DTYPE)

examples/pyseqm_interface/single_node/vars

@@ -0,0 +1,15 @@
+#!/bin/bash
+THIS_PATH=`pwd`
+export PYTHONPATH=""
+MOD_PATH1=$THIS_PATH/../../../src/
+MOD_PATH2=$THIS_PATH/../../../src/sedacs/driver
+PROXYA_FORTRAN_PATH=$THIS_PATH/../../../proxies/fortran/build/
+PROXYA_PYTHON_FILES=$THIS_PATH/../../../proxies/python/
+PROXYA_PYTHON_FOLDER=$THIS_PATH/../../../
+export PROXYA_FORTRAN_PATH=$PROXYA_FORTRAN_PATH
+export PROXYA_PYTHON_PATH=$PROXYA_PYTHON_PATH
+export PYTHONPATH=$PYTHONPATH:$MOD_PATH1:$MOD_PATH2:$PROGRESS_PATH:$BML_PATH:$PROXYA_PYTHON_FOLDER:$PROXYA_PYTHON_FILES:$PROXYA_FORTRAN_PATH
+export PYTHONWARNINGS="ignore:Unverified HTTPS request"
+
+GPU_PATH=$THIS_PATH/../../../src/sedacs/gpu
+export PYTHONPATH=$PYTHONPATH:$GPU_PATH

pyproject.toml

@@ -19,7 +19,7 @@ dependencies = [
     "numpy>=1",
     "torch>=2",
 ]
-requires-python = ">=3.10"
+requires-python = ">=3.9"
 readme = "README.md"
 license = {text = "BSD-3-Clause"}
 classifiers = [
@@ -49,19 +49,19 @@ graph = [
     "metis>=0.2a5"
 ]
 
-[build-system]
-requires = ["scikit-build-core","setuptools", "wheel"]
-build-backend = "scikit_build_core.build"
+#[build-system]
+#requires = ["scikit-build-core","setuptools", "wheel"]
+#build-backend = "scikit_build_core.build"
 
-[tool.scikit-build]
-minimum-version = "0.5"
-cmake.source-dir = "src/sedacs/gpu/nvda/"
-cmake.args = [
-  "-DCMAKE_CXX_COMPILER=nvc++",
-  "-DCMAKE_CXX_FLAGS=-O3 -cuda -gpu=cc89 -acc=gpu -Minfo=accel -fPIC --diag_suppress=bad_macro_redef --diag_suppress=cuda_compile",
-  "-DCMAKE_VERBOSE_MAKEFILE=ON"
-]
-build-dir = "src/sedacs/gpu/nvda/build"
+#[tool.scikit-build]
+#minimum-version = "0.5"
+#cmake.source-dir = "src/sedacs/gpu/nvda/"
+#cmake.args = [
+#  "-DCMAKE_CXX_COMPILER=nvc++",
+#  "-DCMAKE_CXX_FLAGS=-O3 -cuda -gpu=cc89 -acc=gpu -Minfo=accel -fPIC --diag_suppress=bad_macro_redef --diag_suppress=cuda_compile",
+#  "-DCMAKE_VERBOSE_MAKEFILE=ON"
+#]
+#build-dir = "src/sedacs/gpu/nvda/build"
 
 
 [tool.pdm]

src/sedacs/driver/graph_adaptive.py

@@ -412,7 +412,9 @@ def get_singlePoint(sdc,
 
     if rank == 0:
         tic = time.perf_counter()
-        mu0 = get_mu(mu0, full_dVals, full_eVals, sdc.Tel, sy.numel/2)  # oldR
+        #mu0 = get_mu(mu0, full_dVals, full_eVals, sdc.Tel, sy.numel/2)  # oldR
+        mu0 = get_mu(mu0, full_eVals, sdc.Tel, sy.numel/2, dvals = full_dVals)
+
         print("Time mu0 {:>9.4f} (s)".format(time.perf_counter() - tic))
 
     return (EELEC,
@@ -1555,7 +1557,7 @@ def get_adaptiveDM_PYSEQM(sdc,
 
     # Initial chemical potential guess. TODO: This should probably not be 
     # hard-coded.
-    mu0 = -5.5
+    mu0 = -4.7
     if sdc.UHF:
         mu0 = np.array([mu0+0.1, mu0-0.1])
         mu0 = np.array([-1.3, -5.5])
@@ -1682,6 +1684,7 @@ def get_adaptiveDM_PYSEQM(sdc,
                 # are fewer GPUs per node than ranks per nodes.
 
                 if node_rank < num_gpus:
+                    print(mu0)
                     # We want more ranks per node because dm update always
                     # happens on CPU, on node 0, in parallel.
                     (eElec,
@@ -1707,6 +1710,8 @@ def get_adaptiveDM_PYSEQM(sdc,
                                                        graph_for_pairs,
                                                        graph_maskd)
 
+                    print(mu0)
+                    #if gsc > 2: exit(0)
                     gpu_comm.Allreduce(eElec, global_Eelec, op=MPI.SUM)
 
                 else:

src/sedacs/driver/init.py

@@ -122,7 +122,7 @@ def init(args):
 
         # Get hindex (the orbital index for each atom in the system)
         sy.norbs, sy.orbs, hindex, sy.numel, sy.znuc = get_hindex(sdc.orbs, sy.symbols, sy.types)
-        if eng.interface == "PySEQM": sy.numel = int(sy.numel/2)
+        #if eng.interface == "PySEQM": sy.numel = int(sy.numel/2)
         sy.numel -= sdc.charge
         if sdc.UHF:
             sy.nocc_alpha = sy.numel/2. + (sdc.mult-1)/2.

src/sedacs/hamiltonian.py

@@ -96,7 +96,7 @@ def get_hamiltonian(sdc, eng, coords, types, symbols,
         else:
             with torch.no_grad():
                 molSub = get_molecule_pyseqm(sdc, molecule_whole.coordinates[:,partsCoreHaloIndex], symbols, types, device=P_contr.device)[0]
-        M_sub, _, __, ___ = hcore(molSub, doTETCI=False) # off-diagonal h1elec
+        M_sub, _, __, ___, _, _ = hcore(molSub, doTETCI=False) # off-diagonal h1elec
         del _, __, ___
         ham_timing['h1elNonDi'] = time.time() - tic
         #print("t: h1elNonDi {:>7.3f} |".format(time.time() - tic), end=" ")
@@ -213,7 +213,7 @@ def get_hamiltonian(sdc, eng, coords, types, symbols,
             #print("idxi&idxj {:>7.3f} |".format(time.time() - tic), end=" ")
 
             tic = time.time() # compute 2c2e and diagonal h1elec
-            coulInts_test, e1b, e2a, _, _ = TETCI(molecule_whole.const, iii, jjj,
+            coulInts_test, e1b, e2a, _, _, _, _ = TETCI(molecule_whole.const, iii, jjj,
                     molecule_whole.Z[iii], molecule_whole.Z[jjj], x_ij, r_ij, molecule_whole.Z,\
                     molecule_whole.parameters['zeta_s'], molecule_whole.parameters['zeta_p'], molecule_whole.parameters['zeta_d'],\
                     molecule_whole.parameters['s_orb_exp_tail'], molecule_whole.parameters['p_orb_exp_tail'], molecule_whole.parameters['d_orb_exp_tail'],\