From b4257ae2d82bf42495592da4459f24536b0f9ca4 Mon Sep 17 00:00:00 2001 From: Daniel Pozsar Date: Tue, 22 Oct 2024 13:26:44 +0200 Subject: [PATCH] added test.ipynb --- .gitignore | 1 - test.ipynb | 722 +++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 722 insertions(+), 1 deletion(-) create mode 100644 test.ipynb diff --git a/.gitignore b/.gitignore index fb3d223..70a852a 100644 --- a/.gitignore +++ b/.gitignore @@ -21,7 +21,6 @@ docs/source/api/generated # Debug tmp* .coverage -test.ipynb # Mac stuff *.DS_Store* diff --git a/test.ipynb b/test.ipynb new file mode 100644 index 0000000..b817f6c --- /dev/null +++ b/test.ipynb @@ -0,0 +1,722 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "import os\n", + "os.environ[\"OMP_NUM_THREADS\"] = \"4\" # export OMP_NUM_THREADS=4\n", + "os.environ[\"OPENBLAS_NUM_THREADS\"] = \"4\" # export OPENBLAS_NUM_THREADS=4 \n", + "os.environ[\"MKL_NUM_THREADS\"] = \"4\" # export MKL_NUM_THREADS=6\n", + "os.environ[\"VECLIB_MAXIMUM_THREADS\"] = \"4\" # export VECLIB_MAXIMUM_THREADS=4\n", + "os.environ[\"NUMEXPR_NUM_THREADS\"] = \"4\" # export NUMEXPR_NUM_THREADS=6" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "from scipy.linalg import eigh\n", + "#import matplotlib.pyplot as plt\n", + "import warnings\n", + "import sisl\n", + "from grogu.useful import *" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "'0.15.1'" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "np.set_printoptions(linewidth=2*75)\n", + "sisl.__version__" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "# this cell mimicks an input file\n", + "fdf = sisl.get_sile('./Fe3GeTe2/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 sfor this\n", + "ref_xcf_orientations=[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", + "# human readable definition of magnetic entities\n", + "magnetic_entities=[dict(atom=3 ,l=2),\n", + " dict(atom=4 ,l=2),\n", + " dict(atom=5 ,l=2),\n", + " dict(atom=[3,4],)]\n", + "\n", + "# pair information\n", + "pairs=[dict(ai=0,aj=1,Ruc=np.array([0,0,0])), # isotropic should be -82 meV\n", + " dict(ai=0,aj=2,Ruc=np.array([0,0,0])), # 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", + "# Bz sampling\n", + "nk=np.array([10,10,0])" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "ename": "FileNotFoundError", + "evalue": "[Errno 2] No such file or directory: 'Fe3GeTe2/monolayer/soc/lat3_791/Fe3GeTe2.fdf'", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mFileNotFoundError\u001b[0m Traceback (most recent call last)", + "Cell \u001b[0;32mIn[5], line 3\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;66;03m# digestion of the input\u001b[39;00m\n\u001b[1;32m 2\u001b[0m \u001b[38;5;66;03m# read in geometry and hamiltonian\u001b[39;00m\n\u001b[0;32m----> 3\u001b[0m geo \u001b[38;5;241m=\u001b[39m \u001b[43mfdf\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mread_geometry\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 4\u001b[0m dh \u001b[38;5;241m=\u001b[39m fdf\u001b[38;5;241m.\u001b[39mread_hamiltonian()\n\u001b[1;32m 6\u001b[0m uc_in_sc_idx\u001b[38;5;241m=\u001b[39mdh\u001b[38;5;241m.\u001b[39mlattice\u001b[38;5;241m.\u001b[39msc_index([\u001b[38;5;241m0\u001b[39m,\u001b[38;5;241m0\u001b[39m,\u001b[38;5;241m0\u001b[39m])\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/siesta/fdf.py:1440\u001b[0m, in \u001b[0;36mfdfSileSiesta.read_geometry\u001b[0;34m(self, output, *args, **kwargs)\u001b[0m\n\u001b[1;32m 1434\u001b[0m order \u001b[38;5;241m=\u001b[39m _parse_order(\n\u001b[1;32m 1435\u001b[0m kwargs\u001b[38;5;241m.\u001b[39mpop(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124morder\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;28;01mNone\u001b[39;00m),\n\u001b[1;32m 1436\u001b[0m {\u001b[38;5;28;01mTrue\u001b[39;00m: [\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mXV\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mnc\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mTSHS\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mSTRUCT\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mHSX\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mfdf\u001b[39m\u001b[38;5;124m\"\u001b[39m], \u001b[38;5;28;01mFalse\u001b[39;00m: \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mfdf\u001b[39m\u001b[38;5;124m\"\u001b[39m},\n\u001b[1;32m 1437\u001b[0m output,\n\u001b[1;32m 1438\u001b[0m )\n\u001b[1;32m 1439\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m f \u001b[38;5;129;01min\u001b[39;00m order:\n\u001b[0;32m-> 1440\u001b[0m v \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mgetattr\u001b[39;49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43mf\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43m_r_geometry_\u001b[39;49m\u001b[38;5;132;43;01m{\u001b[39;49;00m\u001b[43mf\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mlower\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\u001b[38;5;132;43;01m}\u001b[39;49;00m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43margs\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mkwargs\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 1441\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m v \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m:\n\u001b[1;32m 1442\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m v\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/siesta/fdf.py:1511\u001b[0m, in \u001b[0;36mfdfSileSiesta._r_geometry_fdf\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 1506\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21m_r_geometry_fdf\u001b[39m(\u001b[38;5;28mself\u001b[39m, \u001b[38;5;241m*\u001b[39margs, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs):\n\u001b[1;32m 1507\u001b[0m \u001b[38;5;250m \u001b[39m\u001b[38;5;124;03m\"\"\"Returns Geometry object from the FDF file\u001b[39;00m\n\u001b[1;32m 1508\u001b[0m \n\u001b[1;32m 1509\u001b[0m \u001b[38;5;124;03m NOTE: Interaction range of the Atoms are currently not read.\u001b[39;00m\n\u001b[1;32m 1510\u001b[0m \u001b[38;5;124;03m \"\"\"\u001b[39;00m\n\u001b[0;32m-> 1511\u001b[0m lattice \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mread_lattice\u001b[49m\u001b[43m(\u001b[49m\u001b[43morder\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mfdf\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[1;32m 1513\u001b[0m \u001b[38;5;66;03m# No fractional coordinates\u001b[39;00m\n\u001b[1;32m 1514\u001b[0m is_frac \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;01mFalse\u001b[39;00m\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/siesta/fdf.py:774\u001b[0m, in \u001b[0;36mfdfSileSiesta.read_lattice\u001b[0;34m(self, output, *args, **kwargs)\u001b[0m\n\u001b[1;32m 768\u001b[0m order \u001b[38;5;241m=\u001b[39m _parse_order(\n\u001b[1;32m 769\u001b[0m kwargs\u001b[38;5;241m.\u001b[39mpop(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124morder\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;28;01mNone\u001b[39;00m),\n\u001b[1;32m 770\u001b[0m {\u001b[38;5;28;01mTrue\u001b[39;00m: [\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mXV\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mnc\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mTSHS\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mSTRUCT\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mHSX\u001b[39m\u001b[38;5;124m\"\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mfdf\u001b[39m\u001b[38;5;124m\"\u001b[39m], \u001b[38;5;28;01mFalse\u001b[39;00m: \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mfdf\u001b[39m\u001b[38;5;124m\"\u001b[39m},\n\u001b[1;32m 771\u001b[0m output,\n\u001b[1;32m 772\u001b[0m )\n\u001b[1;32m 773\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m f \u001b[38;5;129;01min\u001b[39;00m order:\n\u001b[0;32m--> 774\u001b[0m v \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mgetattr\u001b[39;49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43mf\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43m_r_lattice_\u001b[39;49m\u001b[38;5;132;43;01m{\u001b[39;49;00m\u001b[43mf\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mlower\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\u001b[38;5;132;43;01m}\u001b[39;49;00m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43margs\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mkwargs\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 775\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m v \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m:\n\u001b[1;32m 776\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m v\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/siesta/fdf.py:781\u001b[0m, in \u001b[0;36mfdfSileSiesta._r_lattice_fdf\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 779\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21m_r_lattice_fdf\u001b[39m(\u001b[38;5;28mself\u001b[39m, \u001b[38;5;241m*\u001b[39margs, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs):\n\u001b[1;32m 780\u001b[0m \u001b[38;5;250m \u001b[39m\u001b[38;5;124;03m\"\"\"Returns `Lattice` object from the FDF file\"\"\"\u001b[39;00m\n\u001b[0;32m--> 781\u001b[0m s \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mget\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mLatticeConstant\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43munit\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mAng\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[1;32m 782\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m s \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m:\n\u001b[1;32m 783\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m SileError(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mCould not find LatticeConstant in file\u001b[39m\u001b[38;5;124m\"\u001b[39m)\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/sile.py:707\u001b[0m, in \u001b[0;36msile_fh_open.._wrapper..pre_open\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 705\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mhasattr\u001b[39m(\u001b[38;5;28mself\u001b[39m, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mfh\u001b[39m\u001b[38;5;124m\"\u001b[39m):\n\u001b[1;32m 706\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m func(\u001b[38;5;28mself\u001b[39m, \u001b[38;5;241m*\u001b[39margs, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs)\n\u001b[0;32m--> 707\u001b[0m \u001b[38;5;28;01mwith\u001b[39;00m \u001b[38;5;28mself\u001b[39m:\n\u001b[1;32m 708\u001b[0m reset(\u001b[38;5;28mself\u001b[39m)\n\u001b[1;32m 709\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m func(\u001b[38;5;28mself\u001b[39m, \u001b[38;5;241m*\u001b[39margs, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs)\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/sile.py:1107\u001b[0m, in \u001b[0;36mSile.__enter__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 1105\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21m__enter__\u001b[39m(\u001b[38;5;28mself\u001b[39m):\n\u001b[1;32m 1106\u001b[0m \u001b[38;5;250m \u001b[39m\u001b[38;5;124;03m\"\"\"Opens the output file and returns it self\"\"\"\u001b[39;00m\n\u001b[0;32m-> 1107\u001b[0m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_open\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 1108\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28mself\u001b[39m\n", + "File \u001b[0;32m~/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.12/site-packages/sisl/io/sile.py:1098\u001b[0m, in \u001b[0;36mSile._open\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 1096\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mfh \u001b[38;5;241m=\u001b[39m gzip\u001b[38;5;241m.\u001b[39mopen(\u001b[38;5;28mstr\u001b[39m(\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mfile), mode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_mode)\n\u001b[1;32m 1097\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[0;32m-> 1098\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mfh \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfile\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopen\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_mode\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 1100\u001b[0m \u001b[38;5;66;03m# the file should restart the file-read (as per instructed)\u001b[39;00m\n\u001b[1;32m 1101\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_line \u001b[38;5;241m=\u001b[39m \u001b[38;5;241m0\u001b[39m\n", + "File \u001b[0;32m/opt/homebrew/Cellar/python@3.12/3.12.7_1/Frameworks/Python.framework/Versions/3.12/lib/python3.12/pathlib.py:1013\u001b[0m, in \u001b[0;36mPath.open\u001b[0;34m(self, mode, buffering, encoding, errors, newline)\u001b[0m\n\u001b[1;32m 1011\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mb\u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;129;01min\u001b[39;00m mode:\n\u001b[1;32m 1012\u001b[0m encoding \u001b[38;5;241m=\u001b[39m io\u001b[38;5;241m.\u001b[39mtext_encoding(encoding)\n\u001b[0;32m-> 1013\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[43mio\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mopen\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mmode\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mbuffering\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mencoding\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43merrors\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mnewline\u001b[49m\u001b[43m)\u001b[49m\n", + "\u001b[0;31mFileNotFoundError\u001b[0m: [Errno 2] No such file or directory: 'Fe3GeTe2/monolayer/soc/lat3_791/Fe3GeTe2.fdf'" + ] + } + ], + "source": [ + "# digestion of the input\n", + "# read in geometry and hamiltonian\n", + "geo = fdf.read_geometry()\n", + "dh = fdf.read_hamiltonian()\n", + "\n", + "uc_in_sc_idx=dh.lattice.sc_index([0,0,0])\n", + "\n", + "# get indecies for orbitals of the magnetic entities\n", + "for i,d in enumerate(magnetic_entities):\n", + " parsed = parse_magnetic_entity(dh,**d)\n", + " magnetic_entities[i]['orbital_indeces'] = parsed\n", + " magnetic_entities[i]['spin_box_indeces'] = blow_up_orbindx(parsed)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# WE WILL NOT NEED THIS!!\n", + "eigfile=sisl.io.siesta.eigSileSiesta('./Fe3GeTe2/Fe3GeTe2.EIG')\n", + "EF=eigfile.read_fermi_level()\n", + "EF" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "NO=dh.no # shorthand for number of orbitals in the unit cell\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 = dh.tocsr(dh.S_idx).toarray().reshape(NO,dh.n_s,NO).transpose(0,2,1).astype('complex128')\n", + "\n", + "# Reorganization of Hamiltonian and overlap matrix elements to SPIN BOX representation\n", + "\n", + "U=np.vstack([np.kron(np.eye(NO,dtype=int),[1,0]),np.kron(np.eye(NO,dtype=int),[0,1])])\n", + "# This is the permutation that transforms ud1ud2 to u12d12\n", + "# That is this transforms FROM SPIN BOX to ORBITAL BOX => U\n", + "\n", + "# the inverse transformation is U.T u12d12 to ud1ud2\n", + "# That is FROM ORBITAL BOX to SPIN BOX => U.T\n", + "\n", + "# this is a test\n", + "# u12d12=np.array([np.kron([1,1],np.arange(NO)),np.kron([0,1],np.ones(NO))],dtype=int).T\n", + "# ud1ud2=np.array([np.kron(np.arange(NO),[1,1]),np.kron(np.ones(NO),[0,1])],dtype=int).T\n", + "# np.allclose(u12d12,U@ud1ud2)\n", + "\n", + "hh,ss = \\\n", + "np.array([U.T@np.block([[h11[:,:,i],h12[:,:,i]],\n", + " [h21[:,:,i],h22[:,:,i]]])@U for i in range(dh.lattice.nsc.prod())]), \\\n", + "np.array([U.T@np.block([[sov[:,:,i] ,sov[:,:,i]*0],\n", + " [sov[:,:,i]*0,sov[:,:,i] ]])@U for i in range(dh.lattice.nsc.prod())])\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())]\n", + "\n", + "XCF = np.array([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", + "# 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(f\"Exchange field has non negligeble scalr part. Largest value is {max_xcfs}\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# collecting information for all proposed rotations, this might generate a lot of unnceessary data\n", + "# information beeing collected: \n", + "# - matrices rotating the direction of the exchange field to the proposed reference direction\n", + "# - rotated xc field \n", + "# - single commutators of the unit cell localized exchange field for the two variations perpendicular to the reference direction\n", + "# - double commutators of the unit cell localized exchange field for the two variations perpendicular to the reference direction\n", + "\n", + "\n", + "# Transformation matrices for rotating the exchange field\n", + "for i,orient in enumerate(ref_xcf_orientations):\n", + "\n", + " # obtain rotated exchange field \n", + " R=RotMa2b(scf_xcf_orientation,orient['o'])\n", + " \n", + " rot_XCF = np.einsum('ij,jklm->iklm',R, XCF)\n", + " rot_H_XCF = sum([np.kron(rot_XCF[i],tau) for i,tau in enumerate([tau_x,tau_y,tau_z])])\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 = hTRS+rot_H_XCF\n", + " \n", + " # Update magnetic entities with commutator data for each \n", + " for mag_ent in magnetic_entities:\n", + " mag_ent['Commutator_Data'] = []\n", + "\n", + " for u in orient['vw']:\n", + " \n", + " Tu = np.kron(np.eye(NO,dtype=int),tau_u(u))\n", + "\n", + " Vu1 = 1j/2*comm(rot_H_XCF_uc,Tu)\n", + " Vu2 = 1/8*comm(comm(Tu,rot_H_XCF_uc),Tu)\n", + "\n", + " for mag_ent in magnetic_entities:\n", + " idx = mag_ent['spin_box_indeces']\n", + " mag_ent['Commutator_Data'].append(dict(Vu1=Vu1[idx,:][:,idx],\n", + " Vu2=Vu2[idx,:][:,idx]))\n", + "\n", + "\n", + "# TODO STORE ABOVE INFO!!!!!" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "k = np.array([0.11441,0.234432,0.0])\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "#This goes inside the for loop for k \n", + "k = np.array([0,0.1,0.0])\n", + "\n", + "\n", + "phases = np.exp(-2j*np.pi*k@dh.lattice.sc_off.T)\n", + "\n", + "np.exp(-1j * np.dot(np.dot(np.dot(dh.rcell, k), dh.cell),dh.lattice.sc_off.T))\n", + "ss.shape,rot_H_XCF.shape,hTRS.shape,phases.shape\n", + "\n", + "HK=np.einsum('abc,a->bc',hTRS+rot_H_XCF,phases)\n", + "SK=np.einsum('abc,a->bc',ss,phases)\n", + "\n", + "vals,vecs=eigh(HK,SK)\n", + "\n", + "#This goes inside the for loop for ze\n", + "ze = -.1+0.000j\n", + "GK=vecs@np.diag(1/(ze-vals))@vecs.conj().T\n", + "\n", + "#phase=np.exp(1j*np.dot(np.dot(k,dh.rcell),offset))\n", + "Gii=GK[magn_ent_idices[ai],:][:,magn_ent_idices[ai]]\n", + "Gjj=GK[magn_ent_idices[ai],:][:,magn_ent_idices[ai]]\n", + "Gij=GK[magn_ent_idices[ai],:][:,magn_ent_idices[ai]]*np.exp(-k@offset*2*np.pi)\n", + "Gji=GK[magn_ent_idices[ai],:][:,magn_ent_idices[ai]]*np.exp( k@offset*2*np.pi)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "np.allclose(np.linalg.inv(ze*SK-HK),GK), abs(np.linalg.inv(ze*SK-HK)-GK).max()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def hsk(dh,k=(0,0,0)):\n", + " '''\n", + " One way to speed up Hk and Sk generation\n", + " '''\n", + " k = np.asarray(k, np.float64) # this two conversion lines\n", + " k.shape = (-1,) # are from the sisl source\n", + "\n", + " # this generates the list of phases\n", + " phases = np.exp(-1j * np.dot(np.dot(np.dot(dh.rcell, k), dh.cell),\n", + " dh.lattice.sc_off.T))\n", + "\n", + " HK11 = np.einsum('abc,c->ab',dh.h11,phases)\n", + " HK12 = np.einsum('abc,c->ab',dh.h12,phases)\n", + " HK21 = np.einsum('abc,c->ab',dh.h21,phases)\n", + " HK22 = np.einsum('abc,c->ab',dh.h22,phases)\n", + "\n", + " SK = np.einsum('abc,c->ab',dh.sov,phases)\n", + "\n", + " return HK11,HK12,HK21,HK22,SK\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + " \n", + "\n", + "rot_xcf_mats=[]\n", + "for i,rot in enumerate(rots):\n", + " R=RotM(*rot) # rotation matrix to go from scf to reference direction\n", + " R@scf_xcf_orientation # exchange field orientation \n", + "\n", + "\n", + "rot_xcf_mats = [RotM(*rot) for rot in rots] \n", + "# directions to ratate around \n", + "dirs_xcf = [R@scf_xcf_orientation for R in rot_xcf_mats]\n", + "dirs_deriv=[]\n", + "for o in dirs_xcf:\n", + " x,y,z = np.eye(3)\n", + " R = RotMa2b(z,o)\n", + " dirs_deriv.append([R@x,R@y])\n", + "\n", + "\n", + "# geting the exchange field part of the reference Hamiltonians for the rotated exchange field\n", + "dh.hxc_rots=[]\n", + "for mtx in rot_xcf_mats:\n", + " rotated_xcf = np.einsum('ij,jklm->iklm',mtx, XCF)\n", + " dh.hxc_rots.append(sum([np.kron(rotated_xcf[i],tau) \n", + " for i,tau in enumerate([tau_x,tau_y,tau_z]) ]))\n", + "\n", + "\n", + "\n", + "cummutators=[]\n", + "\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "i=0\n", + "\n", + "dh.hxc_rots[i][uc_in_sc_idx]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dh.hTRS.shape" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "np.kron(np.random.rand(3,3),tau_x)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dh.hh" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Sanity check for the exchange field orientations\n", + "xcf_test=np.array([[abs(xcfield[i].c).max(),\n", + " abs(xcfield[i].x).max(),\n", + " abs(xcfield[i].y).max(),\n", + " abs(xcfield[i].z).max()] for i in range(dh.lattice.nsc.prod())])\n", + "# Check how the exchange field components behave\n", + "plt.plot(np.linalg.norm(dh.lattice.sc_off,axis=1),xcf_test[:,0],'o',label='scalar')\n", + "plt.plot(np.linalg.norm(dh.lattice.sc_off,axis=1),xcf_test[:,1],'o',label='x')\n", + "plt.plot(np.linalg.norm(dh.lattice.sc_off,axis=1),xcf_test[:,2],'o',label='y')\n", + "plt.plot(np.linalg.norm(dh.lattice.sc_off,axis=1),xcf_test[:,3],'o',label='z')\n", + "plt.yscale('log')\n", + "plt.ylim(1e-22,1e1)\n", + "plt.legend()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "i=41\n", + "DAT=np.vstack([np.hstack([dh.h11[:,:,i]+EF*dh.sov[:,:,i],dh.h12[:,:,i]]),\n", + " np.hstack([dh.h21[:,:,i], dh.h22[:,:,i]]+EF*dh.sov[:,:,i])])\n", + "dh.lattice.sc_off[i]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "for i,v in enumerate(dh.lattice.sc_off):\n", + " if np.allclose(v,np.array([1,0,0])):\n", + " print(i)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.pcolor((U.T@DAT@U)[:6,:6].imag)\n", + "plt.colorbar()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dh.h11 = dh.tocsr(dh.M11r)\n", + "dh.h11 += dh.tocsr(dh.M11i)*1.0j\n", + "dh.h11 = dh.h11.toarray().astype('complex128')\n", + "dh.h22 = dh.tocsr(dh.M22r)\n", + "dh.h22 += dh.tocsr(dh.M22i)*1.0j\n", + "dh.h22 = dh.h22.toarray().astype('complex128')\n", + "\n", + "dh.h12 = dh.tocsr(dh.M12r)\n", + "dh.h12 += dh.tocsr(dh.M12i)*1.0j\n", + "dh.h12 = dh.h12.toarray().astype('complex128')\n", + "dh.h21 = dh.tocsr(dh.M21r)\n", + "dh.h21 += dh.tocsr(dh.M21i)*1.0j\n", + "dh.h21 = dh.h21.toarray().reshape(dh.no,dh.n_s,dh.no).transpose(0,2,1).astype('complex128')\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dh.tocsr(dh.S_idx" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def hsk(dh,k=(0,0,0)):\n", + " '''\n", + " One way to speed up Hk and Sk generation\n", + " '''\n", + " k = np.asarray(k, np.float64) # this two conversion lines\n", + " k.shape = (-1,) # are from the sisl source\n", + "\n", + " # this generates the list of phases\n", + " phases = np.exp(-1j * np.dot(np.dot(np.dot(dh.rcell, k), dh.cell),\n", + " dh.lattice.sc_off.T))\n", + "\n", + " HK11 = np.einsum('abc,c->ab',dh.h11,phases)\n", + " HK12 = np.einsum('abc,c->ab',dh.h12,phases)\n", + " HK21 = np.einsum('abc,c->ab',dh.h21,phases)\n", + " HK22 = np.einsum('abc,c->ab',dh.h22,phases)\n", + "\n", + " SK = np.einsum('abc,c->ab',dh.sov,phases)\n", + "\n", + " return HK11,HK12,HK21,HK22,SK\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "myHk=(lambda x:np.vstack([np.hstack([x[0],x[1]]),\n", + " np.hstack([x[2],x[3]])]))(hsk(dh))\n", + "mySk=(lambda x:np.vstack([np.hstack([x[-1],0*x[-1]]),\n", + " np.hstack([0*x[-1],x[-1]])]))(hsk(dh))\n", + "\n", + "myHkshift=myHk+EF*mySk" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "H0=np.vstack( [np.hstack([dh.h11[:,:,0],dh.h12[:,:,0]]),\n", + " np.hstack([dh.h21[:,:,0],dh.h22[:,:,0]])])\n", + "\n", + "H0shifted=np.vstack( [np.hstack([dh.h11[:,:,0]+EF*dh.sov[:,:,0],dh.h12[:,:,0]]),\n", + " np.hstack([dh.h21[:,:,0], dh.h22[:,:,0]+EF*dh.sov[:,:,0]])])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "(U.T@H0shifted@U)[:3,:3]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.pcolor((U.T@H0shifted@U)[:5,:5])\n", + "plt.colorbar()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "np.linalg.eigvalsh(myHk+EF*mySk)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dh.Hk().toarray()[:3,:3]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.matshow(np.abs(dh.Hk().toarray()-myHk))\n", + "plt.colorbar()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.plot(np.linalg.eigvalsh(dh.Hk().toarray())-np.linalg.eigvalsh(myHk))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "magnetic_entities = [parse_magnetic_entity(dh,atom=[0,1],l=2),parse_magnetic_entity(dh,atom=2,l=2)]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# TODO\n", + "\n", + "- Definition of magnetic entities:\n", + " * Through simple sequence o forbitals in the unit cell\n", + " * Through atom specification\n", + " * Through atom and orbital specification\n", + "- Separation of TR and TRB components of the Hamiltonian, Identification of the exchange field. \n", + "- Definition of commutator expressions, old projection matrix elements\n", + "- Efficient calculation of Green's functions\n", + "- Calculation of energy and momentum resolved derivatives\n", + "- Parallel BZ and serial energy integral" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "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.12.7" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +}