diff --git a/Fe3GeTe2_notebook.pickle b/Fe3GeTe2_notebook.pickle new file mode 100644 index 0000000..922da06 Binary files /dev/null and b/Fe3GeTe2_notebook.pickle differ diff --git a/README.md b/README.md index 6871ccf..058e3c2 100644 --- a/README.md +++ b/README.md @@ -3,7 +3,12 @@ More on the theoretical background can be seen on [arXiv](https://arxiv.org/abs/ # TODO ## Testing -- Run tests on different magnetic materials and compare it to Grogu Matlab --> ran on Jij_for_Marci_6p45ang, but I could not compare data +- Test on lat3_791 +- Test on varcell +- Test on varcell_myrun +- Test on varcell_myrun_2 +- Test scaling +- Test rotations to x-y 30 and 60 ## Developing - Magnetic entities cannot handle orbitals, only shells @@ -12,6 +17,7 @@ More on the theoretical background can be seen on [arXiv](https://arxiv.org/abs/ - Check if exchange field has scalar part - Add more tests!! - io stuff +- logging # Building wheel See detailed documentation on [PYPI](https://packaging.python.org/en/latest/tutorials/packaging-projects/). diff --git a/src/grogupy/__init__.py b/src/grogupy/__init__.py index c055582..ce4a624 100644 --- a/src/grogupy/__init__.py +++ b/src/grogupy/__init__.py @@ -23,6 +23,7 @@ from grogupy.core import * +from grogupy.globals import * from grogupy.io import * from grogupy.magnetism import * from grogupy.utilities import * diff --git a/src/grogupy/globals.py b/src/grogupy/globals.py new file mode 100644 index 0000000..9b76a27 --- /dev/null +++ b/src/grogupy/globals.py @@ -0,0 +1,61 @@ +# Copyright (c) [2024] [] +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in all +# copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE +# SOFTWARE. + +import numpy as np + +# Pauli matrices +TAU_X = np.array([[0, 1], [1, 0]]) +TAU_Y = np.array([[0, -1j], [1j, 0]]) +TAU_Z = np.array([[1, 0], [0, -1]]) +TAU_0 = np.array([[1, 0], [0, 1]]) + + +# list of accepted input parameters +ACCEPTED_INPUTS = [ + "infile", + "outfile", + "scf_xcf_orientation", + "ref_xcf_orientations", + "kset", + "kdirs", + "ebot", + "eset", + "esetp", + "parallel_solver_for_Gk", + "padawan_mode", +] + +DEFAULT_ARGUMENTS = dict( + infile=None, + outfile=None, + scf_xcf_orientation=np.array([0, 0, 1]), + ref_xcf_orientations=[ + dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]), + dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), + dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), + ], + kset=2, + kdirs="xyz", + ebot=None, + eset=42, + esetp=1000, + parallel_solver_for_Gk=False, + padawan_mode=True, +) diff --git a/src/grogupy/grogu.py b/src/grogupy/grogu.py index be77b4c..cc457c8 100644 --- a/src/grogupy/grogu.py +++ b/src/grogupy/grogu.py @@ -72,10 +72,6 @@ def parse_command_line(): type=str, help="Output file name", ) - - # parser.add_argument('--scf-orientation', dest = 'scf_xcf_orientation', default = None, help = 'Output file name') - # parser.add_argument('--ref-orientation', dest = 'ref_xcf_orientations', default = None, help = 'Output file name') - parser.add_argument( "--kset", dest="kset", @@ -157,21 +153,27 @@ def main(simulation_parameters, magnetic_entities, pairs): except: print("numpy version unknown.") - # rename outfile - if not simulation_parameters["outfile"].endswith(".pickle"): - simulation_parameters["outfile"] += ".pickle" - - # if ebot is not given put it 0.1 eV under the smallest energy - 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.") - else: - simulation_parameters["automatic_ebot"] = False + # rename outfile + if not simulation_parameters["outfile"].endswith(".pickle"): + simulation_parameters["outfile"] += ".pickle" + + # if ebot is not given put it 1 eV under the smallest energy + if simulation_parameters["ebot"] is None: + try: + eigfile = simulation_parameters["infile"][:-3] + "EIG" + simulation_parameters["ebot"] = read_siesta_emin(eigfile) - 1 + simulation_parameters["automatic_ebot"] = True + except: + print("Could not determine ebot.") + print("Parameter was not given and .EIG file was not found.") + else: + simulation_parameters["automatic_ebot"] = False + + # add third direction for off-diagonal anisotropy + for orientation in simulation_parameters["ref_xcf_orientations"]: + o1 = orientation["vw"][0] + o2 = orientation["vw"][1] + orientation["vw"].append((o1 + o2) / np.sqrt(2)) # read sile fdf = sisl.get_sile(simulation_parameters["infile"]) @@ -216,7 +218,7 @@ def main(simulation_parameters, magnetic_entities, pairs): ss[i], ss[j] = (s1 + s1d.T.conj()) / 2, (s1d + s1.T.conj()) / 2 # identifying TRS and TRB parts of the Hamiltonian - TAUY = np.kron(np.eye(NO), tau_y) + TAUY = np.kron(np.eye(NO), TAU_Y) hTR = np.array([TAUY @ hh[i].conj() @ TAUY for i in range(dh.lattice.nsc.prod())]) hTRS = (hh + hTR) / 2 hTRB = (hh - hTR) / 2 @@ -294,7 +296,7 @@ def main(simulation_parameters, magnetic_entities, pairs): R = RotMa2b(simulation_parameters["scf_xcf_orientation"], orient["o"]) rot_XCF = np.einsum("ij,jklm->iklm", R, XCF) rot_H_XCF = sum( - [np.kron(rot_XCF[i], tau) for i, tau in enumerate([tau_x, tau_y, tau_z])] + [np.kron(rot_XCF[i], tau) for i, tau in enumerate([TAU_X, TAU_Y, TAU_Z])] ) rot_H_XCF_uc = rot_H_XCF[uc_in_sc_idx] @@ -466,8 +468,8 @@ def main(simulation_parameters, magnetic_entities, pairs): # calculate magnetic parameters for mag_ent in magnetic_entities: - Kxx, Kyy, Kzz, consistency = calculate_anisotropy_tensor(mag_ent) - mag_ent["K"] = np.array([Kxx, Kyy, Kzz]) * sisl.unit_convert("eV", "meV") + K, consistency = calculate_anisotropy_tensor(mag_ent) + mag_ent["K"] = K * sisl.unit_convert("eV", "meV") mag_ent["K_consistency"] = consistency for pair in pairs: @@ -510,17 +512,17 @@ def main(simulation_parameters, magnetic_entities, pairs): if __name__ == "__main__": # loading parameters # it is not clear how to give grogu.fdf path... - command_line_arguments = parse_command_line() - fdf_arguments, magnetic_entities, pairs = read_fdf(command_line_arguments["infile"]) + # command_line_arguments = parse_command_line() + # fdf_arguments, magnetic_entities, pairs = read_fdf(command_line_arguments["infile"]) # right now we do not use any of these input, but it shows # the order of priority when processing arguments default_arguments = False fdf_arguments = False command_line_arguments = False - simulation_parameters = process_input_args( - default_arguments, fdf_arguments, command_line_arguments, ACCEPTED_INPUTS - ) + # simulation_parameters = process_input_args( + # default_arguments, fdf_arguments, command_line_arguments, ACCEPTED_INPUTS + # ) #################################################################################################### # This is the input file for now # @@ -552,7 +554,7 @@ if __name__ == "__main__": dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), ] - simulation_parameters["kset"] = 20 + simulation_parameters["kset"] = 10 simulation_parameters["kdirs"] = "xy" simulation_parameters["ebot"] = None simulation_parameters["eset"] = 600 diff --git a/src/grogupy/io.py b/src/grogupy/io.py index 01bf9f4..5e2aea0 100644 --- a/src/grogupy/io.py +++ b/src/grogupy/io.py @@ -24,41 +24,9 @@ from pickle import dump, load import numpy as np +from globals import ACCEPTED_INPUTS, DEFAULT_ARGUMENTS from sisl.io import fdfSileSiesta -# list of accepted input parameters -ACCEPTED_INPUTS = [ - "infile", - "outfile", - "scf_xcf_orientation", - "ref_xcf_orientations", - "kset", - "kdirs", - "ebot", - "eset", - "esetp", - "parallel_solver_for_Gk", - "padawan_mode", -] - -default_arguments = dict( - infile=None, - outfile=None, - scf_xcf_orientation=np.array([0, 0, 1]), - ref_xcf_orientations=[ - dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]), - dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), - dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), - ], - kset=2, - kdirs="xyz", - ebot=None, - eset=42, - esetp=1000, - parallel_solver_for_Gk=False, - padawan_mode=True, -) - def read_fdf(path): """It reads the simulation parameters, magnetic entities and pairs from the fdf. @@ -397,7 +365,8 @@ def print_atoms_and_pairs(magnetic_entities, pairs): # coordinates and tag print(f"{tag} {xyz[0]} {xyz[1]} {xyz[2]}") print("Consistency check: ", mag_ent["K_consistency"]) - print("Anisotropy diag: ", mag_ent["K"]) + print("K: # Kxx, Kxy, Kxz, Kyx, Kyy, Kyz, Kzx, Kzy, Kzz") + print(mag_ent["K"]) print("") print( diff --git a/src/grogupy/magnetism.py b/src/grogupy/magnetism.py index de58489..0f3d73f 100644 --- a/src/grogupy/magnetism.py +++ b/src/grogupy/magnetism.py @@ -117,8 +117,6 @@ def parse_magnetic_entity(dh, atom=None, l=None, **kwargs): def calculate_anisotropy_tensor(mag_ent): """Calculates the renormalized anisotropy tensor from the energies. - It uses the grogu convention for output. - Args: mag_ent : dict An element from the magnetic entities @@ -126,22 +124,34 @@ def calculate_anisotropy_tensor(mag_ent): Returns: K : np.array_like elements of the anisotropy tensor + consistency_check : float + Absolute value of the difference from the consistency check """ # get the energies energies = mag_ent["energies"] + K = np.zeros((3, 3)) + # calculate the diagonal tensor elements - Kxx = energies[1, 1] - energies[1, 0] - Kyy = energies[0, 1] - energies[0, 0] - Kzz = 0 + K[0, 0] = energies[1, 1] - energies[1, 0] + K[1, 1] = energies[0, 1] - energies[0, 0] + K[2, 2] = 0 + + # calculate the off-diagonal tensor elements + K[0, 1] = (energies[2, 0] + energies[2, 1]) / 2 - energies[2, 2] + K[1, 0] = K[0, 1] + K[0, 2] = (energies[1, 0] + energies[1, 1]) / 2 - energies[1, 2] + K[2, 0] = K[0, 2] + K[1, 2] = (energies[0, 0] + energies[0, 1]) / 2 - energies[0, 2] + K[2, 1] = K[1, 2] # perform consistency check - calculated_diff = Kyy - Kxx + calculated_diff = K[1, 1] - K[0, 0] expected_diff = energies[2, 0] - energies[2, 1] consistency_check = abs(calculated_diff - expected_diff) - return Kxx, Kyy, Kzz, consistency_check + return K, consistency_check def calculate_exchange_tensor(pair): diff --git a/src/grogupy/utilities.py b/src/grogupy/utilities.py index fe63038..e5964d1 100644 --- a/src/grogupy/utilities.py +++ b/src/grogupy/utilities.py @@ -22,15 +22,10 @@ """ import numpy as np +from globals import TAU_0, TAU_X, TAU_Y, TAU_Z from scipy.special import roots_legendre from sisl.io.siesta import eigSileSiesta -# Pauli matrices -tau_x = np.array([[0, 1], [1, 0]]) -tau_y = np.array([[0, -1j], [1j, 0]]) -tau_z = np.array([[1, 0], [0, -1]]) -tau_0 = np.array([[1, 0], [0, 1]]) - def commutator(a, b): """Shorthand for commutator. @@ -189,7 +184,7 @@ def tau_u(u): # u is force to be of unit length u = u / np.linalg.norm(u) - return u[0] * tau_x + u[1] * tau_y + u[2] * tau_z + return u[0] * TAU_X + u[1] * TAU_Y + u[2] * TAU_Z #