et209d
/
grogu
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
Branches Tags
${ item.name }
Create tag ${ searchTerm }
Create branch ${ searchTerm }
from 'cdec373d0e'
${ noResults }
grogu/test.ipynb

810 lines
37 KiB
Raw Blame History

{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"from sys import stdout\n",
"from tqdm import tqdm\n",
"from timeit import default_timer as timer\n",
"\n",
"os.environ[\"OMP_NUM_THREADS\"] = \"1\" # export OMP_NUM_THREADS=4\n",
"os.environ[\"OPENBLAS_NUM_THREADS\"] = \"1\" # export OPENBLAS_NUM_THREADS=4\n",
"os.environ[\"MKL_NUM_THREADS\"] = \"1\" # export MKL_NUM_THREADS=6\n",
"os.environ[\"VECLIB_MAXIMUM_THREADS\"] = \"1\" # export VECLIB_MAXIMUM_THREADS=4\n",
"os.environ[\"NUMEXPR_NUM_THREADS\"] = \"1\" # export NUMEXPR_NUM_THREADS=6\n",
"\n",
"import numpy as np\n",
"import sisl\n",
"from src.grogu_magn.useful import *\n",
"from mpi4py import MPI\n",
"from numpy.linalg import inv\n",
"import warnings\n",
"\n",
"# runtime information\n",
"times = dict()\n",
"times[\"start_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# this cell mimicks an input file\n",
"fdf = sisl.get_sile(\"./Jij_for_Marci_6p45ang/CrBr.fdf\") # ./lat3_791/Fe3GeTe2.fdf\n",
"# this information needs to be given at the input!!\n",
"scf_xcf_orientation = np.array([0, 0, 1]) # z\n",
"# list of reference directions for around which we calculate the derivatives\n",
"# o is the quantization axis, v and w are two axes perpendicular to it\n",
"# at this moment the user has to supply o,v,w on the input.\n",
"# we can have some default for this\n",
"ref_xcf_orientations = [\n",
" dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]),\n",
" dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]),\n",
" dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]),\n",
"]\n",
"\n",
"\"\"\"\n",
"# human readable definition of magnetic entities ./lat3_791/Fe3GeTe2.fdf\n",
"magnetic_entities = [\n",
" dict(atom=3, l=2),\n",
" dict(atom=4, l=2),\n",
" dict(atom=5, l=2),\n",
" dict(\n",
" atom=[3, 4],\n",
" ),\n",
"]\n",
"# pair information ./lat3_791/Fe3GeTe2.fdf\n",
"pairs = [\n",
" dict(ai=0, aj=1, Ruc=np.array([0, 0, 0])), # isotropic should be -82 meV\n",
" dict(\n",
" ai=0, aj=2, Ruc=np.array([0, 0, 0])\n",
" ), # these should all be around -41.9 in the isotropic part\n",
" dict(ai=1, aj=2, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([-1, 0, 0])),\n",
" dict(ai=1, aj=2, Ruc=np.array([-1, 0, 0])),\n",
"] \"\"\"\n",
"\n",
"# human readable definition of magnetic entities ./Jij_for_Marci_6p45ang/CrBr.fdf\n",
"magnetic_entities = [\n",
" dict(atom=0, l=2),\n",
" dict(atom=1, l=2),\n",
" dict(atom=2, l=2),\n",
"]\n",
"# pair information ./Jij_for_Marci_6p45ang/CrBr.fdf\n",
"pairs = [\n",
" dict(ai=0, aj=1, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=1, aj=2, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=0, aj=1, Ruc=np.array([1, 0, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([1, 0, 0])),\n",
" dict(ai=0, aj=1, Ruc=np.array([-1, 0, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([-1, 0, 0])),\n",
" dict(ai=0, aj=1, Ruc=np.array([0, 1, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([0, 1, 0])),\n",
" dict(ai=0, aj=1, Ruc=np.array([0, 1, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([0, 1, 0])),\n",
"]\n",
"\n",
"\n",
"# Brilloun zone sampling and Green function contour integral\n",
"kset = 20\n",
"kdirs = \"xy\"\n",
"ebot = -30\n",
"eset = 100\n",
"esetp = 10000\n",
"\n",
"\n",
"# MPI parameters\n",
"comm = MPI.COMM_WORLD\n",
"size = comm.Get_size()\n",
"rank = comm.Get_rank()\n",
"root_node = 0\n",
"if rank == root_node:\n",
" print(\"Number of nodes in the parallel cluster: \", size)\n",
"\n",
"simulation_parameters = dict(\n",
" path=\"Not yet specified.\",\n",
" scf_xcf_orientation=scf_xcf_orientation,\n",
" ref_xcf_orientations=ref_xcf_orientations,\n",
" kset=kset,\n",
" kdirs=kdirs,\n",
" ebot=ebot,\n",
" eset=eset,\n",
" esetp=esetp,\n",
" parallel_size=size,\n",
")\n",
"\n",
"# digestion of the input\n",
"# read in hamiltonian\n",
"dh = fdf.read_hamiltonian()\n",
"try:\n",
" simulation_parameters[\"geom\"] = fdf.read_geometry()\n",
"except:\n",
" print(\"Error reading geometry.\")\n",
"\n",
"# unit cell index\n",
"uc_in_sc_idx = dh.lattice.sc_index([0, 0, 0])\n",
"\n",
"times[\"setup_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"NO = dh.no # shorthand for number of orbitals in the unit cell\n",
"\n",
"# preprocessing Hamiltonian and overlap matrix elements\n",
"h11 = dh.tocsr(dh.M11r)\n",
"h11 += dh.tocsr(dh.M11i) * 1.0j\n",
"h11 = h11.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype(\"complex128\")\n",
"\n",
"h22 = dh.tocsr(dh.M22r)\n",
"h22 += dh.tocsr(dh.M22i) * 1.0j\n",
"h22 = h22.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype(\"complex128\")\n",
"\n",
"h12 = dh.tocsr(dh.M12r)\n",
"h12 += dh.tocsr(dh.M12i) * 1.0j\n",
"h12 = h12.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype(\"complex128\")\n",
"\n",
"h21 = dh.tocsr(dh.M21r)\n",
"h21 += dh.tocsr(dh.M21i) * 1.0j\n",
"h21 = h21.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype(\"complex128\")\n",
"\n",
"sov = (\n",
" dh.tocsr(dh.S_idx)\n",
" .toarray()\n",
" .reshape(NO, dh.n_s, NO)\n",
" .transpose(0, 2, 1)\n",
" .astype(\"complex128\")\n",
")\n",
"\n",
"\n",
"# Reorganization of Hamiltonian and overlap matrix elements to SPIN BOX representation\n",
"U = np.vstack(\n",
" [np.kron(np.eye(NO, dtype=int), [1, 0]), np.kron(np.eye(NO, dtype=int), [0, 1])]\n",
")\n",
"# This is the permutation that transforms ud1ud2 to u12d12\n",
"# That is this transforms FROM SPIN BOX to ORBITAL BOX => U\n",
"# the inverse transformation is U.T u12d12 to ud1ud2\n",
"# That is FROM ORBITAL BOX to SPIN BOX => U.T\n",
"\n",
"# From now on everything is in SPIN BOX!!\n",
"hh, ss = np.array(\n",
" [\n",
" U.T @ np.block([[h11[:, :, i], h12[:, :, i]], [h21[:, :, i], h22[:, :, i]]]) @ U\n",
" for i in range(dh.lattice.nsc.prod())\n",
" ]\n",
"), np.array(\n",
" [\n",
" U.T\n",
" @ np.block([[sov[:, :, i], sov[:, :, i] * 0], [sov[:, :, i] * 0, sov[:, :, i]]])\n",
" @ U\n",
" for i in range(dh.lattice.nsc.prod())\n",
" ]\n",
")\n",
"\n",
"\n",
"# symmetrizing Hamiltonian and overlap matrix to make them hermitian\n",
"for i in range(dh.lattice.sc_off.shape[0]):\n",
" j = dh.lattice.sc_index(-dh.lattice.sc_off[i])\n",
" h1, h1d = hh[i], hh[j]\n",
" hh[i], hh[j] = (h1 + h1d.T.conj()) / 2, (h1d + h1.T.conj()) / 2\n",
" s1, s1d = ss[i], ss[j]\n",
" ss[i], ss[j] = (s1 + s1d.T.conj()) / 2, (s1d + s1.T.conj()) / 2\n",
"\n",
"# identifying TRS and TRB parts of the Hamiltonian\n",
"TAUY = np.kron(np.eye(NO), tau_y)\n",
"hTR = np.array([TAUY @ hh[i].conj() @ TAUY for i in range(dh.lattice.nsc.prod())])\n",
"hTRS = (hh + hTR) / 2\n",
"hTRB = (hh - hTR) / 2\n",
"\n",
"# extracting the exchange field\n",
"traced = [spin_tracer(hTRB[i]) for i in range(dh.lattice.nsc.prod())] # equation 77\n",
"XCF = np.array(\n",
" [\n",
" np.array([f[\"x\"] for f in traced]),\n",
" np.array([f[\"y\"] for f in traced]),\n",
" np.array([f[\"z\"] for f in traced]),\n",
" ]\n",
") # equation 77\n",
"\n",
"# Check if exchange field has scalar part\n",
"max_xcfs = abs(np.array(np.array([f[\"c\"] for f in traced]))).max()\n",
"if max_xcfs > 1e-12:\n",
" warnings.warn(\n",
" f\"Exchange field has non negligible scalar part. Largest value is {max_xcfs}\"\n",
" )\n",
"\n",
"times[\"H_and_XCF_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# for every site we have to store 3 Greens function (and the associated _tmp-s) in the 3 reference directions\n",
"for i, mag_ent in enumerate(magnetic_entities):\n",
" parsed = parse_magnetic_entity(dh, **mag_ent) # parse orbital indexes\n",
" magnetic_entities[i][\"orbital_indeces\"] = parsed\n",
" # calculate spin box indexes\n",
" magnetic_entities[i][\"spin_box_indeces\"] = blow_up_orbindx(parsed)\n",
" # calculate size for Greens function generation\n",
" spin_box_shape = len(mag_ent[\"spin_box_indeces\"])\n",
"\n",
" mag_ent[\"energies\"] = [] # we will store the second order energy derivations here\n",
"\n",
" mag_ent[\"Gii\"] = [] # Greens function\n",
" mag_ent[\"Gii_tmp\"] = [] # Greens function for parallelization\n",
" # These will be the perturbed potentials from eq. 100\n",
" mag_ent[\"Vu1\"] = [list([]) for _ in range(len(ref_xcf_orientations))]\n",
" mag_ent[\"Vu2\"] = [list([]) for _ in range(len(ref_xcf_orientations))]\n",
" for i in ref_xcf_orientations:\n",
" # Greens functions for every quantization axis\n",
" mag_ent[\"Gii\"].append(\n",
" np.zeros((eset, spin_box_shape, spin_box_shape), dtype=\"complex128\")\n",
" )\n",
" mag_ent[\"Gii_tmp\"].append(\n",
" np.zeros((eset, spin_box_shape, spin_box_shape), dtype=\"complex128\")\n",
" )\n",
"\n",
"# for every site we have to store 2x3 Greens function (and the associated _tmp-s)\n",
"# in the 3 reference directions, because G_ij and G_ji are both needed\n",
"for pair in pairs:\n",
" # calculate size for Greens function generation\n",
" spin_box_shape_i = len(magnetic_entities[pair[\"ai\"]][\"spin_box_indeces\"])\n",
" spin_box_shape_j = len(magnetic_entities[pair[\"aj\"]][\"spin_box_indeces\"])\n",
"\n",
" pair[\"energies\"] = [] # we will store the second order energy derivations here\n",
"\n",
" pair[\"Gij\"] = [] # Greens function\n",
" pair[\"Gji\"] = []\n",
" pair[\"Gij_tmp\"] = [] # Greens function for parallelization\n",
" pair[\"Gji_tmp\"] = []\n",
" for i in ref_xcf_orientations:\n",
" # Greens functions for every quantization axis\n",
" pair[\"Gij\"].append(\n",
" np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype=\"complex128\")\n",
" )\n",
" pair[\"Gij_tmp\"].append(\n",
" np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype=\"complex128\")\n",
" )\n",
" pair[\"Gji\"].append(\n",
" np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype=\"complex128\")\n",
" )\n",
" pair[\"Gji_tmp\"].append(\n",
" np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype=\"complex128\")\n",
" )\n",
"\n",
"times[\"site_and_pair_dictionaries_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"kset = make_kset(dirs=kdirs, NUMK=kset) # generate k space sampling\n",
"wkset = np.ones(len(kset)) / len(kset) # generate weights for k points\n",
"kpcs = np.array_split(kset, size) # split the k points based on MPI size\n",
"kpcs[root_node] = tqdm(kpcs[root_node], desc=\"k loop\", file=stdout)\n",
"\n",
"times[\"k_set_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# this will contain the three hamiltonians in the reference directions needed to calculate the energy variations upon rotation\n",
"hamiltonians = []\n",
"\n",
"# iterate over the reference directions (quantization axes)\n",
"for i, orient in enumerate(ref_xcf_orientations):\n",
" # obtain rotated exchange field\n",
" R = RotMa2b(scf_xcf_orientation, orient[\"o\"])\n",
" rot_XCF = np.einsum(\"ij,jklm->iklm\", R, XCF)\n",
" rot_H_XCF = sum(\n",
" [np.kron(rot_XCF[i], tau) for i, tau in enumerate([tau_x, tau_y, tau_z])]\n",
" )\n",
" rot_H_XCF_uc = rot_H_XCF[uc_in_sc_idx]\n",
"\n",
" # obtain total Hamiltonian with the rotated exchange field\n",
" rot_H = (\n",
" hTRS + rot_H_XCF\n",
" ) # equation 76 #######################################################################################\n",
"\n",
" hamiltonians.append(\n",
" dict(orient=orient[\"o\"], H=rot_H)\n",
" ) # store orientation and rotated Hamiltonian\n",
"\n",
" # these are the infinitezimal rotations (for now) perpendicular to the quantization axis\n",
" for u in orient[\"vw\"]:\n",
" Tu = np.kron(np.eye(NO, dtype=int), tau_u(u)) # section 2.H\n",
"\n",
" Vu1 = 1j / 2 * commutator(rot_H_XCF_uc, Tu) # equation 100\n",
" Vu2 = 1 / 8 * commutator(commutator(Tu, rot_H_XCF_uc), Tu) # equation 100\n",
"\n",
" for mag_ent in magnetic_entities:\n",
" # fill up the perturbed potentials (for now) based on the on-site projections\n",
" mag_ent[\"Vu1\"][i].append(\n",
" Vu1[:, mag_ent[\"spin_box_indeces\"]][mag_ent[\"spin_box_indeces\"], :]\n",
" )\n",
" mag_ent[\"Vu2\"][i].append(\n",
" Vu2[:, mag_ent[\"spin_box_indeces\"]][mag_ent[\"spin_box_indeces\"], :]\n",
" )\n",
"\n",
"times[\"reference_rotations_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"if rank == root_node:\n",
" print(\"Number of magnetic entities being calculated: \", len(magnetic_entities))\n",
" print(\n",
" \"We have to calculate the Greens function for three reference direction and we are going to calculate 15 energy integrals per site.\"\n",
" )\n",
" print(f\"The shape of the Hamiltonian and the Greens function is {NO}x{NO}.\")\n",
"comm.Barrier()\n",
"# ----------------------------------------------------------------------\n",
"\n",
"# make energy contour\n",
"# we are working in eV now !\n",
"# and sisil shifts E_F to 0 !\n",
"cont = make_contour(emin=ebot, enum=eset, p=esetp)\n",
"eran = cont.ze\n",
"\n",
"# ----------------------------------------------------------------------\n",
"# sampling the integrand on the contour and the BZ\n",
"for k in kpcs[rank]:\n",
" wk = wkset[rank] # weight of k point in BZ integral\n",
" # iterate over reference directions\n",
" for i, hamiltonian_orientation in enumerate(hamiltonians):\n",
" # calculate Greens function\n",
" H = hamiltonian_orientation[\"H\"]\n",
" HK, SK = hsk(H, ss, dh.sc_off, k)\n",
" Gk = inv(SK * eran.reshape(eset, 1, 1) - HK)\n",
"\n",
" # solve Greens function sequentially for the energies, because of memory bound\n",
" # Gk = np.zeros(shape=(eset, HK.shape[0], HK.shape[1]), dtype=\"complex128\")\n",
" # for j in range(eset):\n",
" # Gk[j] = inv(SK * eran[j] - HK)\n",
"\n",
" # store the Greens function slice of the magnetic entities (for now) based on the on-site projections\n",
" for mag_ent in magnetic_entities:\n",
" mag_ent[\"Gii_tmp\"][i] += (\n",
" Gk[:, mag_ent[\"spin_box_indeces\"]][..., mag_ent[\"spin_box_indeces\"]]\n",
" * wk\n",
" )\n",
"\n",
" for pair in pairs:\n",
" # add phase shift based on the cell difference\n",
" phase = np.exp(1j * 2 * np.pi * k @ pair[\"Ruc\"].T)\n",
"\n",
" # get the pair orbital sizes from the magnetic entities\n",
" ai = magnetic_entities[pair[\"ai\"]][\"spin_box_indeces\"]\n",
" aj = magnetic_entities[pair[\"aj\"]][\"spin_box_indeces\"]\n",
"\n",
" # store the Greens function slice of the magnetic entities (for now) based on the on-site projections\n",
" pair[\"Gij_tmp\"][i] += Gk[:, ai][..., aj] * phase * wk\n",
" pair[\"Gji_tmp\"][i] += Gk[:, aj][..., ai] * phase * wk\n",
"\n",
"# summ reduce partial results of mpi nodes\n",
"for i in range(len(hamiltonians)):\n",
" for mag_ent in magnetic_entities:\n",
" comm.Reduce(mag_ent[\"Gii_tmp\"][i], mag_ent[\"Gii\"][i], root=root_node)\n",
"\n",
" for pair in pairs:\n",
" comm.Reduce(pair[\"Gij_tmp\"][i], pair[\"Gij\"][i], root=root_node)\n",
" comm.Reduce(pair[\"Gji_tmp\"][i], pair[\"Gji\"][i], root=root_node)\n",
"\n",
"times[\"green_function_inversion_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"if rank == root_node:\n",
" # iterate over the magnetic entities\n",
" for tracker, mag_ent in enumerate(magnetic_entities):\n",
" # iterate over the quantization axes\n",
" for i, Gii in enumerate(mag_ent[\"Gii\"]):\n",
" storage = []\n",
" # iterate over the first and second order local perturbations\n",
" for Vu1, Vu2 in zip(mag_ent[\"Vu1\"][i], mag_ent[\"Vu2\"][i]):\n",
" # The Szunyogh-Lichtenstein formula\n",
" traced = np.trace((Vu2 @ Gii + 0.5 * Gii @ Vu1 @ Gii), axis1=1, axis2=2)\n",
" # evaluation of the contour integral\n",
" storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we)))\n",
"\n",
" # fill up the magnetic entities dictionary with the energies\n",
" magnetic_entities[tracker][\"energies\"].append(storage)\n",
"\n",
" # iterate over the pairs\n",
" for tracker, pair in enumerate(pairs):\n",
" # iterate over the quantization axes\n",
" for i, (Gij, Gji) in enumerate(zip(pair[\"Gij\"], pair[\"Gji\"])):\n",
" site_i = magnetic_entities[pair[\"ai\"]]\n",
" site_j = magnetic_entities[pair[\"aj\"]]\n",
"\n",
" storage = []\n",
" # iterate over the first order local perturbations in all possible orientations for the two sites\n",
" for Vui in site_i[\"Vu1\"][i]:\n",
" for Vuj in site_j[\"Vu1\"][i]:\n",
" # The Szunyogh-Lichtenstein formula\n",
" traced = np.trace((Vui @ Gij @ Vuj @ Gji), axis1=1, axis2=2)\n",
" # evaluation of the contour integral\n",
" storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we)))\n",
" # fill up the pairs dictionary with the energies\n",
" pairs[tracker][\"energies\"].append(storage)\n",
"\n",
" times[\"end_time\"] = timer()\n",
" print_output(simulation_parameters, magnetic_entities, pairs, dh, times)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"##################################################################### GROGU OUTPUT #############################################################################\n",
"================================================================================================================================================================\n",
"Input file: \n",
"Not yet specified.\n",
"Number of nodes in the parallel cluster: 1\n",
"================================================================================================================================================================\n",
"Cell [Ang]: \n",
"[[ 6.47 0. 0. ]\n",
" [-3.235 5.60318436 0. ]\n",
" [ 0. 0. 29.999449 ]]\n",
"================================================================================================================================================================\n",
"DFT axis: \n",
"[0 0 1]\n",
"Quantization axis and perpendicular rotation directions:\n",
"[1 0 0] --» [array([0, 1, 0]), array([0, 0, 1])]\n",
"[0 1 0] --» [array([1, 0, 0]), array([0, 0, 1])]\n",
"[0 0 1] --» [array([1, 0, 0]), array([0, 1, 0])]\n",
"================================================================================================================================================================\n",
"number of k points: 4\n",
"k point directions: xy\n",
"================================================================================================================================================================\n",
"Parameters for the contour integral:\n",
"Ebot: -30\n",
"Eset: 100\n",
"Esetp: 10000\n",
"================================================================================================================================================================\n",
"Atomic informations: \n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[atom index]Element(orbitals) x [Ang] y [Ang] z [Ang] Sx Sy Sz Q Lx Ly Lz Jx Jy Jz\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[0]Cr(2) 0.025224941687312104 -0.014371429816522993 -0.020792425897793128\n",
"\n",
"[1]Cr(2) 3.2508367285399955 1.8478445793345957 -0.02082986959108256\n",
"\n",
"[2]Br(2) 2.1053343349045215 0.10721337338160676 -1.471910834765145\n",
"\n",
"================================================================================================================================================================\n",
"Exchange [meV]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Magnetic entity1 Magnetic entity2 [i j k] d [Ang]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[0]Cr(2) [1]Cr(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: 15.91078131425419\n",
"DMI: [-3.16312278e-05 1.94021443e-05 -1.65972633e-04]\n",
"Symmetric-anisotropy: [-0.03896181 -0.01392426 -0.07639227 0.05293223 -0.08712209]\n",
"Energies for debugging: \n",
"array([[ 1.58968571e-02, 8.70904601e-05, 8.71537225e-05,\n",
" 1.59636674e-02],\n",
" [ 1.58999738e-02, -5.29322330e-05, -5.28934287e-05,\n",
" 1.58718195e-02],\n",
" [ 1.59589510e-02, 7.62263019e-05, 7.65582472e-05,\n",
" 1.58707914e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([0.01587182, 0.01595895, 0.01589686])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"0.015871819503030533 0.015870791381430095\n",
"\n",
"[0]Cr(2) [2]Br(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.14694676346924165\n",
"DMI: [ 8.01967166e-05 8.29086493e-04 -9.27235617e-05]\n",
"Symmetric-anisotropy: [ 2.66339812e-05 -4.77282294e-06 -6.60708923e-06 8.03002228e-04\n",
" -2.98623982e-06]\n",
"Energies for debugging: \n",
"array([[ 1.46941991e-04, 8.31829565e-08, -7.72104768e-08,\n",
" 1.46924902e-04],\n",
" [ 1.46982775e-04, -8.03002228e-07, 8.55170759e-07,\n",
" 1.46973397e-04],\n",
" [ 1.46942514e-04, -8.61164724e-08, 9.93306509e-08,\n",
" 1.46901824e-04]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([0.00014697, 0.00014694, 0.00014694])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"0.00014697339745047312 0.00014690182371456808\n",
"\n",
"[1]Cr(2) [2]Br(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.146966596490555\n",
"DMI: [ 7.83469069e-04 -3.25183011e-04 9.29353203e-05]\n",
"Symmetric-anisotropy: [-1.73027923e-05 1.26068464e-05 -2.08606748e-05 -2.97793770e-04\n",
" 3.14336799e-07]\n",
"Energies for debugging: \n",
"array([[ 1.46979203e-04, 7.83154732e-07, -7.83783406e-07,\n",
" 1.46971292e-04],\n",
" [ 1.46912608e-04, 2.97793770e-07, -3.52572253e-07,\n",
" 1.46949294e-04],\n",
" [ 1.46905570e-04, 1.13795995e-07, -7.20746454e-08,\n",
" 1.46914411e-04]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([0.00014695, 0.00014691, 0.00014698])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"0.00014694929369829223 0.00014691441062702917\n",
"\n",
"[0]Cr(2) [1]Cr(2) [1 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.09111814992141927\n",
"DMI: [0.03377305 0.0340363 0.00964244]\n",
"Symmetric-anisotropy: [-0.00291957 0.0031014 -0.0021492 0.03416682 0.0051237 ]\n",
"Energies for debugging: \n",
"array([[ 9.42195473e-05, 2.86493440e-05, -3.88967487e-05,\n",
" 9.09363256e-05],\n",
" [ 9.41093226e-05, -3.41668178e-05, 3.39057913e-05,\n",
" 8.81985768e-05],\n",
" [ 9.14933866e-05, 1.17916416e-05, -7.49323868e-06,\n",
" 8.86473282e-05]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([8.81985768e-05, 9.14933866e-05, 9.42195473e-05])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"8.819857682132204e-05 8.8647328224418e-05\n",
"\n",
"[0]Cr(2) [2]Br(2) [1 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.00023346314999400802\n",
"DMI: [-0.0002294 -0.0002349 -0.00020879]\n",
"Symmetric-anisotropy: [-1.12556482e-05 1.43988058e-05 -1.55475448e-06 -2.34959571e-04\n",
" 6.19663196e-06]\n",
"Energies for debugging: \n",
"array([[ 2.47861956e-07, -2.35597865e-07, 2.23204601e-07,\n",
" 2.30319992e-07],\n",
" [ 2.10261125e-07, 2.34959571e-07, -2.34837442e-07,\n",
" 2.22207502e-07],\n",
" [ 2.12323581e-07, -2.07238961e-07, 2.10348470e-07,\n",
" 2.07858965e-07]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([2.22207502e-07, 2.12323581e-07, 2.47861956e-07])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"2.2220750178916023e-07 2.078589648170586e-07\n",
"\n",
"[0]Cr(2) [1]Cr(2) [-1 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.09111258017019819\n",
"DMI: [-0.03378783 -0.03403353 -0.0096277 ]\n",
"Symmetric-anisotropy: [-0.00291741 0.00310155 -0.00214665 -0.03390027 0.00512057]\n",
"Energies for debugging: \n",
"array([[ 9.42141309e-05, -3.89083994e-05, 2.86672518e-05,\n",
" 9.09284386e-05],\n",
" [ 9.41042657e-05, 3.39002650e-05, -3.41667971e-05,\n",
" 8.81951711e-05],\n",
" [ 9.14629700e-05, -7.48105019e-06, 1.17743519e-05,\n",
" 8.86174723e-05]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([8.81951711e-05, 9.14629700e-05, 9.42141309e-05])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"8.819517110814147e-05 8.861747229207156e-05\n",
"\n",
"[0]Cr(2) [2]Br(2) [-1 0 0] d [Ang] Not yet.\n",
"Isotropic: 0.0002516119318884681\n",
"DMI: [0.0002271 0.00022722 0.00023639]\n",
"Symmetric-anisotropy: [ 8.12786238e-06 -1.33724299e-05 2.28353801e-06 2.31211900e-04\n",
" -2.16488405e-06]\n",
"Energies for debugging: \n",
"array([[ 2.38239502e-07, 2.29269073e-07, -2.24939305e-07,\n",
" 2.56856499e-07],\n",
" [ 2.77408158e-07, -2.31211900e-07, 2.23219979e-07,\n",
" 2.59739794e-07],\n",
" [ 2.73816864e-07, 2.34111415e-07, -2.38678491e-07,\n",
" 2.71887108e-07]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([2.59739794e-07, 2.73816864e-07, 2.38239502e-07])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"2.5973979426862243e-07 2.718871078150382e-07\n",
"\n",
"[0]Cr(2) [1]Cr(2) [0 1 0] d [Ang] Not yet.\n",
"Isotropic: -0.13657014924460614\n",
"DMI: [ 0.00026859 -0.00029338 -0.0001818 ]\n",
"Symmetric-anisotropy: [-0.00104734 0.00154555 -0.00046486 -0.00072287 0.00087507]\n",
"Energies for debugging: \n",
"array([[-1.35024601e-04, -6.06475771e-07, -1.14365755e-06,\n",
" -1.37068359e-04],\n",
" [-1.34947686e-04, 7.22871593e-07, 1.36112752e-07,\n",
" -1.37617487e-04],\n",
" [-1.36691877e-04, 2.83065496e-07, 6.46663000e-07,\n",
" -1.37228934e-04]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.00013762, -0.00013669, -0.00013502])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.00013761748740593208 -0.00013722893403347796\n",
"\n",
"[0]Cr(2) [2]Br(2) [0 1 0] d [Ang] Not yet.\n",
"Isotropic: -0.00010028627777987032\n",
"DMI: [-8.54009417e-05 -1.21477411e-04 -8.32475011e-05]\n",
"Symmetric-anisotropy: [ 1.89621993e-05 -7.40845210e-06 6.91965022e-07 -1.13803644e-04\n",
" -7.00808812e-06]\n",
"Energies for debugging: \n",
"array([[-1.07694730e-07, -7.83928536e-08, 9.24090298e-08,\n",
" -1.11840025e-07],\n",
" [-8.13196124e-08, 1.13803644e-07, -1.29151178e-07,\n",
" -8.13240785e-08],\n",
" [-1.76542958e-07, -8.39394661e-08, 8.25555360e-08,\n",
" -1.76140161e-07]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-8.13240785e-08, -1.76542958e-07, -1.07694730e-07])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-8.132407846090182e-08 -1.7614016112304844e-07\n",
"\n",
"[0]Cr(2) [1]Cr(2) [0 1 0] d [Ang] Not yet.\n",
"Isotropic: -0.13657014924460614\n",
"DMI: [ 0.00026859 -0.00029338 -0.0001818 ]\n",
"Symmetric-anisotropy: [-0.00104734 0.00154555 -0.00046486 -0.00072287 0.00087507]\n",
"Energies for debugging: \n",
"array([[-1.35024601e-04, -6.06475771e-07, -1.14365755e-06,\n",
" -1.37068359e-04],\n",
" [-1.34947686e-04, 7.22871593e-07, 1.36112752e-07,\n",
" -1.37617487e-04],\n",
" [-1.36691877e-04, 2.83065496e-07, 6.46663000e-07,\n",
" -1.37228934e-04]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.00013762, -0.00013669, -0.00013502])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.00013761748740593208 -0.00013722893403347796\n",
"\n",
"[0]Cr(2) [2]Br(2) [0 1 0] d [Ang] Not yet.\n",
"Isotropic: -0.00010028627777987032\n",
"DMI: [-8.54009417e-05 -1.21477411e-04 -8.32475011e-05]\n",
"Symmetric-anisotropy: [ 1.89621993e-05 -7.40845210e-06 6.91965022e-07 -1.13803644e-04\n",
" -7.00808812e-06]\n",
"Energies for debugging: \n",
"array([[-1.07694730e-07, -7.83928536e-08, 9.24090298e-08,\n",
" -1.11840025e-07],\n",
" [-8.13196124e-08, 1.13803644e-07, -1.29151178e-07,\n",
" -8.13240785e-08],\n",
" [-1.76542958e-07, -8.39394661e-08, 8.25555360e-08,\n",
" -1.76140161e-07]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-8.13240785e-08, -1.76542958e-07, -1.07694730e-07])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-8.132407846090182e-08 -1.7614016112304844e-07\n",
"\n",
"================================================================================================================================================================\n",
"Runtime information: \n",
"Total runtime: 16.641370917 s\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Initial setup: 0.07859033299999996 s\n",
"Hamiltonian conversion and XC field extraction: 0.683 s\n",
"Pair and site datastructure creatrions: 0.012 s\n",
"k set cration and distribution: 0.023 s\n",
"Rotating XC potential: 0.243 s\n",
"Greens function inversion: 15.524 s\n",
"Calculate energies and magnetic components: 0.078 s"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"========================================\n",
" \n",
"Atom Angstrom\n",
"# Label, x y z Sx Sy Sz #Q Lx Ly Lz Jx Jy Jz\n",
"--------------------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Te1 1.8955 1.0943 13.1698 -0.0000 0.0000 -0.1543 # 5.9345 -0.0000 0.0000 -0.0537 -0.0000 0.0000 -0.2080 \n",
"Te2 1.8955 1.0943 7.4002 0.0000 -0.0000 -0.1543 # 5.9345 0.0000 -0.0000 -0.0537 0.0000 -0.0000 -0.2080 \n",
"Ge3 -0.0000 2.1887 10.2850 0.0000 0.0000 -0.1605 # 3.1927 -0.0000 0.0000 0.0012 0.0000 0.0000 -0.1593 \n",
"Fe4 -0.0000 0.0000 11.6576 0.0001 -0.0001 2.0466 # 8.3044 0.0000 -0.0000 0.1606 0.0001 -0.0001 2.2072 \n",
"Fe5 -0.0000 0.0000 8.9124 -0.0001 0.0001 2.0466 # 8.3044 -0.0000 0.0000 0.1606 -0.0001 0.0001 2.2072 \n",
"Fe6 1.8955 1.0944 10.2850 0.0000 0.0000 1.5824 # 8.3296 -0.0000 -0.0000 0.0520 -0.0000 0.0000 1.6344 \n",
"==================================================================================================================================\n",
" \n",
"Exchange meV\n",
"--------------------------------------------------------------------------------\n",
"# at1 at2 i j k # d (Ang)\n",
"--------------------------------------------------------------------------------\n",
"Fe4 Fe5 0 0 0 # 2.7452\n",
"Isotropic -82.0854\n",
"DMI 0.12557 -0.00082199 6.9668e-08\n",
"Symmetric-anisotropy -0.60237 -0.83842 -0.00032278 -1.2166e-05 -3.3923e-05\n",
"--------------------------------------------------------------------------------\n",
"Fe4 Fe6 0 0 0 # 2.5835\n",
"Isotropic -41.9627\n",
"DMI 1.1205 -1.9532 0.0018386\n",
"Symmetric-anisotropy 0.26007 -0.00013243 0.12977 -0.069979 -0.042066\n",
"--------------------------------------------------------------------------------\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import sisl.viz\n",
"\n",
"dh.geometry.tile(2, 1).plot(axes=\"xy\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"\n",
"coords = dh.xyz[-3:]\n",
"\n",
"\n",
"plt.figure(figsize=(15, 5))\n",
"plt.subplot(131)\n",
"plt.scatter(coords[:, 0], coords[:, 2], color=[\"r\", \"g\", \"b\"])\n",
"plt.xlabel(\"x\")\n",
"plt.ylabel(\"z\")\n",
"plt.subplot(132)\n",
"plt.scatter(coords[:, 1], coords[:, 2], color=[\"r\", \"g\", \"b\"])\n",
"plt.xlabel(\"y\")\n",
"plt.ylabel(\"z\")\n",
"plt.subplot(133)\n",
"plt.scatter(coords[:, 0], coords[:, 1], color=[\"r\", \"g\", \"b\"])\n",
"plt.xlabel(\"x\")\n",
"plt.ylabel(\"y\")\n",
"print(\"xyz[-3:]: red, green, blue\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": ".venv",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Reference in new issue View Git Blame Copy Permalink