added onsite projection, padawan mode, updated printing, sane simulation parameters, reducing number of threads

test.py

@@ -18,6 +18,16 @@
 # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
 # SOFTWARE.
 
+
+import os
+
+os.environ["OMP_NUM_THREADS"] = "1"  # export OMP_NUM_THREADS=1
+os.environ["OPENBLAS_NUM_THREADS"] = "1"  # export OPENBLAS_NUM_THREADS=1
+os.environ["MKL_NUM_THREADS"] = "1"  # export MKL_NUM_THREADS=1
+os.environ["VECLIB_MAXIMUM_THREADS"] = "1"  # export VECLIB_MAXIMUM_THREADS=1
+os.environ["NUMEXPR_NUM_THREADS"] = "1"  # export NUMEXPR_NUM_THREADS=1
+
+
 from timeit import default_timer as timer
 
 # runtime information
@@ -38,7 +48,7 @@ from src.grogu_magn import *
 ######################################################################
 ######################################################################
 
-path = (
+infile = (
     "/Users/danielpozsar/Downloads/nojij/Fe3GeTe2/monolayer/soc/lat3_791/Fe3GeTe2.fdf"
 )
 outfile = "./Fe3GeTe2_notebook"
@@ -60,6 +70,15 @@ pairs = [
     dict(ai=1, aj=2, Ruc=np.array([-3, 0, 0])),
 ]
 simulation_parameters = default_args
+simulation_parameters["infile"] = infile
+simulation_parameters["outfile"] = outfile
+simulation_parameters["kset"] = 20
+simulation_parameters["kdirs"] = "xy"
+simulation_parameters["eset"] = 600
+simulation_parameters["esetp"] = 10000
+
+fdf = sisl.io.fdfSileSiesta("input.fdf")
+fdf.get("XCF_Rotation")
 
 ######################################################################
 ######################################################################
@@ -72,14 +91,17 @@ size = comm.Get_size()
 rank = comm.Get_rank()
 root_node = 0
 
+# include parallel size in simulation parameters
+simulation_parameters["parallel_size"] = size
+
 # check versions for debugging
 if rank == root_node:
     try:
-        print(sisl.__version__)
+        print("sisl version: ", sisl.__version__)
     except:
         print("sisl version unknown.")
     try:
-        print(np.__version__)
+        print("numpy version: ", np.__version__)
     except:
         print("numpy version unknown.")
 
@@ -92,6 +114,7 @@ if simulation_parameters["ebot"] is None:
     try:
         eigfile = simulation_parameters["infile"][:-3] + "EIG"
         simulation_parameters["ebot"] = read_siesta_emin(eigfile) - 0.1
+        simulation_parameters["automatic_ebot"] = True
     except:
         print("Could not determine ebot.")
         print("Parameter was not given and .EIG file was not found.")
@@ -109,7 +132,16 @@ simulation_parameters["cell"] = fdf.read_geometry().cell
 uc_in_sc_idx = dh.lattice.sc_index([0, 0, 0])
 
 if rank == root_node:
-    print_parameters(simulation_parameters)
+    print("\n\n\n\n\n")
+    print(
+        "#################################################################### JOB INFORMATION ###########################################################################"
+    )
+    print_job_description(simulation_parameters)
+    print(
+        "################################################################################################################################################################"
+    )
+    print("\n\n\n\n\n")
+
     times["setup_time"] = timer()
     print(f"Setup done. Elapsed time: {times['setup_time']} s")
     print(
@@ -186,7 +218,7 @@ wkset = np.ones(len(kset)) / len(kset)
 kpcs = np.array_split(kset, size)
 
 # use progress bar if available
-if tqdm_imported:
+if rank == root_node and tqdm_imported:
     kpcs[root_node] = tqdm(kpcs[root_node], desc="k loop")
 
 if rank == root_node:
@@ -236,8 +268,8 @@ for i, orient in enumerate(simulation_parameters["ref_xcf_orientations"]):
         for mag_ent in magnetic_entities:
             idx = mag_ent["spin_box_indices"]
             # fill up the perturbed potentials (for now) based on the on-site projections
-            mag_ent["Vu1"][i].append(onsite_projection(Vu1, idx))
+            mag_ent["Vu1"][i].append(onsite_projection(Vu1, idx, idx))
-            mag_ent["Vu2"][i].append(onsite_projection(Vu2, idx))
+            mag_ent["Vu2"][i].append(onsite_projection(Vu2, idx, idx))
 
 if rank == root_node:
     times["reference_rotations_time"] = timer()
@@ -252,13 +284,17 @@ if rank == root_node:
 # requirements of the Greens function calculations
 if rank == root_node:
     print("Starting matrix inversions.")
+    if simulation_parameters["padawan_mode"]:
+        print("Padawan mode: ")
         print(f"Total number of k points: {kset.shape[0]}")
         print(f"Number of energy samples per k point: {simulation_parameters['eset']}")
         print(f"Total number of directions: {len(hamiltonians)}")
         print(
             f"Total number of matrix inversions: {kset.shape[0] * len(hamiltonians) * simulation_parameters['eset']}"
         )
-    print(f"The shape of the Hamiltonian and the Greens function is {NO}x{NO}={NO*NO}")
+        print(
+            f"The shape of the Hamiltonian and the Greens function is {NO}x{NO}={NO*NO}"
+        )
         # https://stackoverflow.com/questions/70746660/how-to-predict-memory-requirement-for-np-linalg-inv
         # memory is O(64 n**2) for complex matrices
         memory_size = getsizeof(hamiltonians[0]["H"].base) / 1024
@@ -311,7 +347,7 @@ for k in kpcs[rank]:
         # store the Greens function slice of the magnetic entities
         for mag_ent in magnetic_entities:
             idx = mag_ent["spin_box_indices"]
-            mag_ent["Gii_tmp"][i] += onsite_projection(Gk, idx) * wk
+            mag_ent["Gii_tmp"][i] += onsite_projection(Gk, idx, idx) * wk
 
         for pair in pairs:
             # add phase shift based on the cell difference
@@ -411,14 +447,20 @@ if rank == root_node:
         pair["J"] = J * sisl.unit_convert("eV", "meV")
 
     times["end_time"] = timer()
+    print("\n\n\n\n\n")
     print(
         "##################################################################### GROGU OUTPUT #############################################################################"
     )
 
     print_parameters(simulation_parameters)
     print_atoms_and_pairs(magnetic_entities, pairs)
+    print(
+        "################################################################################################################################################################"
+    )
     print_runtime_information(times)
+    print("")
 
+    # remove unwanted stuff before saving
     pairs, magnetic_entities = remove_clutter_for_save(pairs, magnetic_entities)
     # create output dictionary with all the relevant data
     results = dict(
@@ -428,4 +470,7 @@ if rank == root_node:
         runtime=times,
     )
 
+    # save results
     save_pickle(simulation_parameters["outfile"], results)
+
+    print("\n\n\n\n\n")