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.
grogu/test.ipynb

1319 lines
83 KiB

{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"[Mac:23435] shmem: mmap: an error occurred while determining whether or not /var/folders/yh/dx7xl94n3g52ts3td8qcxjcc0000gn/T//ompi.Mac.501/jf.0/46137344/sm_segment.Mac.501.2c00000.0 could be created.\n"
]
},
{
"data": {
"text/plain": [
"'0.14.3'"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"import os\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()\n",
"########################\n",
"# it works if data is in downloads folder\n",
"########################\n",
"sisl.__version__"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"================================================================================================================================================================\n",
"Input file: \n",
"Not yet specified.\n",
"Number of nodes in the parallel cluster: 1\n",
"================================================================================================================================================================\n",
"Cell [Ang]: \n",
"[[ 3.79100000e+00 0.00000000e+00 0.00000000e+00]\n",
" [-1.89550000e+00 3.28310231e+00 0.00000000e+00]\n",
" [ 1.25954923e-15 2.18160327e-15 2.05700000e+01]]\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",
"Parameters for the contour integral:\n",
"Number of k points: 10\n",
"k point directions: xy\n",
"Ebot: -15\n",
"Eset: 50\n",
"Esetp: 1000\n",
"================================================================================================================================================================\n"
]
}
],
"source": [
"# this cell mimicks an input file\n",
"fdf = sisl.get_sile(\n",
" \"/Users/danielpozsar/Downloads/nojij/Fe3GeTe2/monolayer/soc/lat3_791/Fe3GeTe2.fdf\"\n",
") # ./Jij_for_Marci_6p45ang/CrBr.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(atom=[3, 4],),\n",
"]\n",
"# pair information ./lat3_791/Fe3GeTe2.fdf\n",
"pairs = [\n",
" # isotropic should be -82 meV\n",
" dict(ai=0, aj=1, Ruc=np.array([0, 0, 0])),\n",
" # dict(ai=1, aj=0, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=1, aj=2, Ruc=np.array([0, 0, 0])),\n",
" # dict(ai=2, aj=1, Ruc=np.array([0, 0, 0])),\n",
" dict(ai=0, aj=2, Ruc=np.array([0, 0, 0])),\n",
" # dict(ai=2, aj=0, Ruc=np.array([0, 0, 0])),\n",
" # these should all be around -41.9 in the isotropic part\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=2, Ruc=np.array([-1, 0, 0])),\n",
" # dict(ai=1, aj=2, Ruc=np.array([-1, 0, 0])),\n",
"]\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 = 10\n",
"kdirs = \"xy\"\n",
"ebot = -15\n",
"eset = 50\n",
"esetp = 1000\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",
"\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",
"simulation_parameters[\"cell\"] = fdf.read_geometry().cell\n",
"\n",
"# unit cell index\n",
"uc_in_sc_idx = dh.lattice.sc_index([0, 0, 0])\n",
"\n",
"print_parameters(simulation_parameters)\n",
"times[\"setup_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Users/danielpozsar/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.9/site-packages/matplotlib/cbook.py:1762: ComplexWarning: Casting complex values to real discards the imaginary part\n",
" return math.isfinite(val)\n",
"/Users/danielpozsar/Documents/oktatás/elte/phd/grogu_project/.venv/lib/python3.9/site-packages/matplotlib/cbook.py:1398: ComplexWarning: Casting complex values to real discards the imaginary part\n",
" return np.asarray(x, float)\n"
]
},
{
"data": {
"text/plain": [
"-12.806739"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGiCAYAAAAba+fDAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA3T0lEQVR4nO3deXhU9aH/8c9km5BtQvYEEkgAQfZFjHErKrJorSi1bq1g3YWqoK3SX0XBalx6tddeqn1uRbRWqfa6VNuqiCxawhaIC0sEDCRAFkhIJguZJDPn9wdmdCRAApk5M5P363nmeZgzZ5JPjieZj9/zPedYDMMwBAAAECRCzA4AAADQnSg3AAAgqFBuAABAUKHcAACAoEK5AQAAQYVyAwAAggrlBgAABBXKDQAACCqUGwAAEFQoNwAAIKh4tdzk5+dr/Pjxio2NVUpKiqZNm6bi4mKPdZqbmzVr1iwlJiYqJiZG06dPV2VlpTdjAQCAIObVcrNq1SrNmjVLa9eu1bJly9Ta2qpJkyapsbHRvc6cOXP07rvv6o033tCqVau0f/9+XXnlld6MBQAAgpjFlzfOPHDggFJSUrRq1Sqdf/75qqurU3Jysl599VX9+Mc/liRt375dp59+ugoKCnTWWWf5KhoAAAgSYb78ZnV1dZKkhIQESVJhYaFaW1s1ceJE9zpDhgxRVlbWMcuNw+GQw+FwP3e5XKqpqVFiYqIsFouXfwIAANAdDMNQfX29MjIyFBLSvQeSfFZuXC6X7rnnHp1zzjkaPny4JKmiokIRERGKj4/3WDc1NVUVFRUdfp38/HwtWLDA23EBAIAPlJWVqW/fvt36NX1WbmbNmqUvv/xSn3766Sl9nXnz5mnu3Lnu53V1dcrKylJZWZni4uJONSYAdAuXy1De48vV6HDqrTvP1qDUWLMjAX7FbrcrMzNTsbHd/7vhk3Ize/Zsvffee1q9erVHO0tLS1NLS4tqa2s9Rm8qKyuVlpbW4deyWq2yWq1HLY+Li6PcAPAbpdVNOiyrIqNCNDInXeGhXHkD6Ig3ppR49bfNMAzNnj1bb731lj7++GNlZ2d7vD5u3DiFh4dr+fLl7mXFxcUqLS1VXl6eN6MBgFdtLbdLkgalxlBsAB/z6sjNrFmz9Oqrr+qdd95RbGysex6NzWZTr169ZLPZdNNNN2nu3LlKSEhQXFycfvGLXygvL48zpQAEtG3flJvT0xlRBnzNq+XmueeekyRNmDDBY/mLL76omTNnSpKeeeYZhYSEaPr06XI4HJo8ebL++Mc/ejMWAHiVo82p9SU1kig3gBl8ep0bb7Db7bLZbKqrq2PODQBTtTpd+tcX5frdh8UqqzksSfrH7HM0sm+8ucEAP+TNz2+fXucGAIJNq9Ol9z7frw+3VOrTHQdV72iTJKXEWvXA1CEUG8AElBsAOAkul6F3P9+vZ5Z9pd3VTe7lidERmnF2f918XraiIvgTC5iB3zwAOAm3v1KoD7ceuclvYnSErj+rny4ckqKRfWwKCeFq6YCZKDcA0EWtTpc+3l4lSbr7okG69fwcRVv5cwr4C34bAaCL9tceVpvLUERYiO6+aBAjNYCf4cpSANBFe76ZY5OVEEWxAfwQ5QYAumhPzZFy0z8xyuQkADpCuQGALiqtbpQkZSVEm5wEQEcoNwDQRe2HpfoxcgP4JcoNAHRR6TeHpbIoN4BfotwAQBcYhvHtyE0C5QbwR5QbAOiCA/UOHW51KsQi9e1NuQH8EeUGALqg/UypjPheigjjTyjgj/jNBIAuYDIx4P8oNwDQBZwGDvg/yg0AdEH7YSlGbgD/RbkBgC7gTCnA/1FuAKAL9nxzWKpfIoelAH9FuQGATrI3t+pQU6skLuAH+DPKDQB0Uuk3h6SSYiIUYw0zOQ2AY6HcAEAntc+3yWK+DeDXKDcA0El7aphvAwQCyg0AdEJVfbPe2rRPEqeBA/6Og8YAcAJlNU366QvrtKe6ScmxVv14XF+zIwE4DsoNABzH3kNNuur5AlXYm5WZ0Et/veksbpgJ+DnKDQAcQ5vTpXuWFqnC3qxBKTF65eZcpcZFmh0LwAkw5wYAjmHRil3auOeQYqxhWjxzPMUGCBCUGwDoQOGeQ3r24x2SpEemDVMmp38DAYNyAwDfU1HXrLte2yyny9DlozN0xRgmEAOBhHIDAN9R3eDQ9X9eq321h9U/MUqPTBtudiQAXUS5AYBv1B1u1c9eWK9dBxqVbovUKzfnKi4y3OxYALqIcgMAkgzD0F2vbdbWcruSYiL015tzOeUbCFCUGwCQtHRDmVZ9dUDWsBC9/PNc5STHmB0JwEmi3ADo8fYeatJv39sqSfrl5MEamhFnciIAp4JyA6BHc7kM/ervn6uxxakz+vXWjedkmx0JwCmi3ADo0f5euFdrdlUrMjxET101SqEhFrMjAThFlBsAPVaDo01PflAsSZp78WnKToo2ORGA7kC5AdBjPbdypw42ONQ/MUozz+ZwFBAsvFpuVq9ercsuu0wZGRmyWCx6++23PV6fOXOmLBaLx2PKlCnejAQAko5MIv7fT0okSfMuOV0RYfy/HhAsvPrb3NjYqFGjRmnRokXHXGfKlCkqLy93P1577TVvRgIASdKT7xerpc2l3OwETRqaanYcAN0ozJtffOrUqZo6depx17FarUpLS/NmDADwsOtAg/7x2X5ZLNKDPxwqi4VJxEAwMX0cduXKlUpJSdHgwYN1xx13qLq6+rjrOxwO2e12jwcAdMUbG/dKki4cnKLhfWwmpwHQ3UwtN1OmTNHLL7+s5cuX64knntCqVas0depUOZ3OY74nPz9fNpvN/cjMzPRhYgCBrs3p0pubjpSbq87gbt9AMPLqYakTueaaa9z/HjFihEaOHKkBAwZo5cqVuuiiizp8z7x58zR37lz3c7vdTsEB0Gmf7DioqnqHEqIjdOEQ5toAwcj0w1LflZOTo6SkJO3cufOY61itVsXFxXk8AKCz3igskyRdPjqDM6SAIOVXv9l79+5VdXW10tPTzY4CIAjVNLZo2dZKSdJV4xjxBYKVVw9LNTQ0eIzClJSUqKioSAkJCUpISNCCBQs0ffp0paWladeuXfrVr36lgQMHavLkyd6MBaCH+r/CvWp1GhreJ46bYwJBzKvlZuPGjbrgggvcz9vnysyYMUPPPfecPv/8c7300kuqra1VRkaGJk2apEceeURWq9WbsQAEuQZHmz7dcVBV9c1qbnWqyu7Qyq8OaGdVgyRGbYBg59VyM2HCBBmGcczXP/jgA29+ewBBrL65VfXNbTIkHW5xqrSmUSUHm1Swq1qrdxxQS5vrqPeEhVg0eViafnIG5QYIZqaeLQUAXVVlb9Z/L9+hpRvK5HQd+3+e+idGaUhanKIiQhVtDdOZ2Qn6weBkxUWG+zAtADNQbgD4tS/31emVtXvU0uZSi9Ol5duqdLj1yLWwIkJDZLFI4aEhykyIUnZSlE5Pi9OkYWk6LTWGKw8DPRTlBoBfW/DuFm3Yfchj2ZiseD0wZYhycxJNSgXAn1FuAPitwy1OFZXVSpLmTDxN0dZQDUiJ0YTTkhmVAXBMlBsAfmtT6SG1Og1l2CJ110UDKTQAOsWvLuIHAN+17usjN9LNzUmk2ADoNMoNAL+1tqRGkpSbnWByEgCBhHIDwC81t3473+ZMyg2ALqDcAPBLRWW1amlzKTnWquykaLPjAAgglBsAfmnd198ekmK+DYCuoNwA8EvrSr6dTAwAXUG5AeB3Wtpc2lR65MJ9ZzHfBkAXUW4A+J2Ne2rU3OpSYnSEBqbEmB0HQICh3ADwK+9/Wa7b/lIoSTpnYBLzbQB0GVcoBmCaspomrSyuUnFlvVrbDFU3OvTRtipJ0ujMeP36ktNNTgggEFFuAPhUS5tLf9tQqiVrdmvXgcYO17ntBzm6b9JghYcyuAyg6yg3AHziQL1DH26t0B9X7NK+2sOSpNAQi8b1663x/XurV3iowkJDNL5/b43rxyRiACePcgPAqz7YUqE/fLxDX+6zu5elxFo1+8KBmjamj+Iiw01MByAYUW4AeE1tU4t+8epmtThdkqThfeI0bXQfXZ/bT70iQk1OByBYUW4AeM07RfvV4nRpSFqs/nJTrpJjrWZHAtADMFsPgNe8vrFMknTN+EyKDQCfodwA8Iot++u0Zb9dEaEhunx0H7PjAOhBKDcAvOKNjXslSRcPTVXv6AiT0wDoSSg3ALqdo82pt4v2SZKuOqOvyWkA9DSUGwDd7t9fVKi2qVVpcZE6b1Cy2XEA9DCcLQWg27hchhb/p0RPvL9d0pFRm9AQ7g0FwLcoNwC6RXOrU3e8UqgVxQckSZOHpeqOCQNMTgWgJ6LcAOgWS9eXakXxAUWEhWj+D4fq+tws7ugNwBSUGwCnzDAM/e2bs6PmTR2in57Vz+REAHoyJhQDOGVb9tu1rfzINW2uGMM1bQCYi3ID4JS1X4l40rBUxUdxTRsA5qLcADglza1Ovb35yDVtfnJGpslpAIByA+AUfbClQvbmNvWJ76VzBiaZHQcAmFAMoOsOtzhVuOeQSmua9Je1eyRJ08dxTRsA/oFyA6DTGh1t+svaPfrzJ1/rYEOLe3mIRbpqHLdZAOAfKDcAjuuDLRVa/dUB7TrQoC377apvbpMkpcVFamhGnDJ799J5g5KVmRBlclIAOIJyA+CYKu3Nuv2VQhnGt8tykqI164KBunx0hsJCmbYHwP9QbgAcU8nBRhmGlBRj1a8vGaKBKTEalmFjbg0Av+bV/+1avXq1LrvsMmVkZMhisejtt9/2eN0wDM2fP1/p6enq1auXJk6cqB07dngzEoAu2HfosCRpcFqMrhzbVyP7xlNsAPg9r5abxsZGjRo1SosWLerw9SeffFLPPvusnn/+ea1bt07R0dGaPHmympubvRkLQCft/abc9I1nPg2AwOHVw1JTp07V1KlTO3zNMAz9/ve/129+8xtdfvnlkqSXX35Zqampevvtt3XNNdd4MxqATthX2yRJ6tO7l8lJAKDzTJsNWFJSooqKCk2cONG9zGazKTc3VwUFBcd8n8PhkN1u93gA8I59tUdGbvrEU24ABA7Tyk1FRYUkKTU11WN5amqq+7WO5Ofny2azuR+ZmVzuHfAW92EpRm4ABJCAO49z3rx5qqurcz/KysrMjgQEJZfLUHntkflvHJYCEEhMKzdpaWmSpMrKSo/llZWV7tc6YrVaFRcX5/EA0P0ONDjU4nQpNMSitLhIs+MAQKeZVm6ys7OVlpam5cuXu5fZ7XatW7dOeXl5ZsUC8I32Q1JpcZFcrA9AQPHq2VINDQ3auXOn+3lJSYmKioqUkJCgrKws3XPPPfrtb3+rQYMGKTs7Ww8++KAyMjI0bdo0b8YC0Al7D3GmFIDA5NVys3HjRl1wwQXu53PnzpUkzZgxQ0uWLNGvfvUrNTY26tZbb1Vtba3OPfdcvf/++4qMZAgcMFv7mVJ9OVMKQIDxarmZMGGCjO/elOZ7LBaLFi5cqIULF3ozBoCT0H51YkZuAAQaDqQD6BDXuAEQqCg3ADr07TVuuPUCgMBCuQFwFMMwOCwFIGBRbgAc5VBTqw63OiVJ6TYm+AMILJQbAEdpH7VJjrUqMjzU5DQA0DWUGwBHab/GDfeUAhCIKDcAjsKZUgACGeUGwFH2MpkYQACj3AA4ivuwFCM3AAKQV69QDCAw1Da1aE91kz7fV6e/F+7VZ2W1kqS+CVzjBkDgodwAPZRhGHqnaL9+92Gx+zBUu7AQi6aOSFdeTqJJ6QDg5FFugB6ktqlFew8dVmlNk174tESFew65X0uJtap/UrQmDU3VtDF9lBRjNTEpAJw8yg0QxFqdLv35kxJt2F2jL/fVqare4fF6r/BQzb5woGac3V8xVv4cAAgO/DUDgthbm/fpife3eyxLirGqT+9eGp4Rp9kXDlS6jUnDAIIL5QYIYu8U7ZMkXTmmj64/K0tD0uIUzQgNgCDHXzkgSFXam7VmV7Ukac7FpymTM58A9BBc5wYIUu9+tl+GIY3r15tiA6BHodwAQeqdov2SpMtHZ5icBAB8i3IDBKFdBxr0xb46hYZYdOmIdLPjAIBPUW6AIPSPb0ZtzhuUpESuVwOgh2FCMRAEDMPQ3kOHtan0kDbtOaR/fHak3Ewb3cfkZADge5QbIIBtLj2k//3ka60vOaSDDZ4X6EuIjtDFQ1NNSgYA5qHcAAGorKZJT35QrHe/GaGRpPBQi4Zm2DQ2K15js3rr7AGJXNMGQI/EXz4gQLhchv6z66BeWbtHH22rktNlyGKRpo/tq2vGZ2p4H5siw0PNjgkApqPcAAFgy/46/frNL/TZ3jr3snMHJmneJUM0LMNmYjIA8D+UG8CPNbc69cxHX+nPn5TI6TIUHRGq6eP66vrcfhqcFmt2PADwS5QbwE+1OV2a/eomfbStSpJ06Yh0PXTZUKXERZqcDAD8G+UG8EOGYWj+P7boo21VsoaF6Nlrx2jysDSzYwFAQKDcAH6mzenSohW79Oq6Ulks0n9fQ7EBgK6g3AB+wOUy9MnOg/pH0X4t316p2qZWSdLDlw3TlOEUGwDoCsoNYKL65la993m5Xvi0RDurGtzLe0eF6/YfDNCMs/ubFw4AAhTlBvAxp8vQ3wvL9M8vKlSw66BanYYkKcYapulj++iSEeka16+3wkK59RsAnAzKDeBDza1O3b10sz7YUulelpMUretys3T1+EzFRoabmA4AggPlBvCR2qYW3fzSRm3cc0gRoSGafeFAXTIiXQNTYsyOBgBBhXIDeFl53WG9vmGvXltfqgp7s2Ijw/S/N5yhs3ISzY4GAEGJcgN40csFu/XwP7bIdWRajfrE99LimeO5ujAAeBHlBvCSRkebnvqgWC5DOrN/gq7LzdKU4Wnc3BIAvIxyA3jJW5v3qb65TdlJ0Vp661kKCbGYHQkAegTONQW8wDAMvVywW5L0s7P6UWwAwIdMLzcPP/ywLBaLx2PIkCFmxwJOScHX1fqqskFREaH68Rl9zY4DAD2KXxyWGjZsmD766CP387Awv4gFnLSX1uyWJF05to/iuHYNAPiUX7SIsLAwpaV17v45DodDDofD/dxut3srFtBpza1ObdpzSKU1Taqqd2jZ1iMX6ZuR19/cYADQA/lFudmxY4cyMjIUGRmpvLw85efnKysrq8N18/PztWDBAh8nBDwdamxRUVmtNpfVakNJjQpLD6mlzeWxztkDEjUolVO+AcDXLIZhGGYG+Pe//62GhgYNHjxY5eXlWrBggfbt26cvv/xSsbFHfzB0NHKTmZmpuro6xcXF+TI6eqgPt1Tojr9uktPl+auTGmfVsAybkmIilBIbqavHZyozIcqklADg3+x2u2w2m1c+v00fuZk6dar73yNHjlRubq769eun119/XTfddNNR61utVlmtVl9GBDys3nFATpeh1DirzhmQpDFZ8cobkKQBydGyWDgrCgDMZnq5+b74+Hiddtpp2rlzp9lRgA5V1B0ZOfzFhYP007P6mZwGAPB9pp8K/n0NDQ3atWuX0tPTzY4CdKjS3ixJSouLNDkJAKAjppeb++67T6tWrdLu3bu1Zs0aXXHFFQoNDdW1115rdjSgQxXt5cZGuQEAf2T6Yam9e/fq2muvVXV1tZKTk3Xuuedq7dq1Sk5ONjsacJRWp0sHG44clkpl5AYA/JLp5Wbp0qVmRwA6rareIcOQwkMtSoyOMDsOAKADph+WAgJJRd2RQ1IpsZHcLwoA/BTlBuiCSubbAIDfo9wAXdA+csOZUgDgvyg3QBe0j9wwmRgA/BflBuiC8vaRGxtXyQYAf0W5AbqggpEbAPB7lBugC9oPS6XbepmcBABwLJQboJMMw2BCMQAEAMoN0El1h1vlaHNJklLimHMDAP6KcgN0Uvt8m95R4YoMDzU5DQDgWCg3QCe1nynFZGIA8G+UG6CTKuu4OjEABALKDdBJFe4zpSg3AODPKDdAJ3F1YgAIDJQboJM4DRwAAgPlBuikCrtDkpTKYSkA8GuUG6CTKuoOS2LkBgD8HeUG6ITmVqcONbVKotwAgL8LMzsA4M8cbU4VldZqc1mtJMkaFqL4qHBzQwEAjotyAxzDocYWXfWnAu2sanAv658YLYvFYmIqAMCJUG6ADrS0uXT7K4XaWdUgW69wjc2K18CUGF0xpq/Z0QAAJ0C5Ab7HMAz95u0vtK6kRjHWML1+W54Gp8WaHQsA0EmUG/R4DY42ldce1t5Dh7VxT43+s7NaRWW1CrFIf7huDMUGAAIM5QY9UpvTpX98tl/PrdylHd+ZU9POYpEW/GiYLhicYkI6AMCpoNygxyiuqNdne2u1o7JeH2ypVGlNk/u12MgwZdh6aXgfm3JzEpSXk6jMhCgT0wIAThblBj3CH1fu1JPvF3ssS4iO0M3nZeu6M7MUHxVhUjIAQHej3CDovVO0z11scrMTdHp6nIZmxOmHI9MVFcGvAAAEG/6yI6ht2F2jX77xuSTppnOz9eAPh5qcCADgbdx+AUHr/S8rdNOSDWpxujR5WKp+fcnpZkcCAPgAIzcIOo2ONv32n1v12voySdIZ/Xrr91ePUWgIVxYGgJ6AcoOA5XIZsje3qrapVYeaWlRa06QPt1ZqxfYqNbU4ZbFIt50/QHMvPk0RYQxSAkBPQblBQFj11QHl/2ubGlva5HQaamp1qu5wqwyj4/X7JUYp/8oROntAkm+DAgBMR7lBQHh+5S5tr6jv8LXoiFDFR0UoITpCZw9M1CXD0zWyr40bXAJAD0W5gd9ramnTxj01kqQ//Wyc0m2RigwPVXxUuGy9wmUNCzU5IQDAn1Bu4PfWldSo1WmoT3wvTRqayogMAOC4mGUJv/fpjoOSpHMHJlFsAAAnRLmB32svN+edxuRgAMCJUW7g16rszSqurJfFIp3DmU8AgE7wi3KzaNEi9e/fX5GRkcrNzdX69evNjgQ/8ck3ozbDM2zqHc3NLQEAJ2Z6ufnb3/6muXPn6qGHHtKmTZs0atQoTZ48WVVVVWZHgx/4dOc3820GMWoDAOgc08vN008/rVtuuUU33nijhg4dqueff15RUVFavHix2dFgMsMw3OXmvIGUGwBA55hablpaWlRYWKiJEye6l4WEhGjixIkqKCjo8D0Oh0N2u93jgeC0cc8hHah3KDI8ROP69zY7DgAgQJhabg4ePCin06nU1FSP5ampqaqoqOjwPfn5+bLZbO5HZmamL6LCx9aX1OjnSzZIki4cksKF+gAAnWb6Yamumjdvnurq6tyPsrIysyOhm73/ZYV+9sI61Te3aXz/3sq/YqTZkQAAAcTUKxQnJSUpNDRUlZWVHssrKyuVlpbW4XusVqusVqsv4sHHquqb9eg/t+mdov2SpImnp+p/rhujyHBGbQAAnWdquYmIiNC4ceO0fPlyTZs2TZLkcrm0fPlyzZ4928xo8JED9Q5t2F2jdV9X683N+1Tf3CaLRfr5OdmaN3WIwkIDbnARAGAy0+8tNXfuXM2YMUNnnHGGzjzzTP3+979XY2OjbrzxRrOjoZu1tLm0o6pe28rrVbjnkNaVVOvrA40e64zoY9OjVwzXyL7x5oQEAAQ808vN1VdfrQMHDmj+/PmqqKjQ6NGj9f777x81yRiBq6q+WY+8t03//qJcbS7jqNeHpMXqzOwEnT0gURcPTVNoCPePAgCcPIthGEd/2gQQu90um82muro6xcXFmR2nxzMMQ9WNLWp0tKnV6VLhnkN67F/bVXe4VZIUFxmmIelxGtnHptycRI3v31vxUVx5GAB6Gm9+fps+coPA19zq1O8+KNby7VXaV3tYLW2uo9YZ3idOj10xQiP62LizNwDAqyg3OCW7Dzbqzr9u0tZyz4spRkWEKiIsRDHWMP30rH66+dxsJgcDAHyCcoOTtuqrA5r1101qcLQpITpCCy8fplF945UaF6mIMIoMAMAclBuclKKyWt3+l0IdbnXqzP4JevbaMUqzRZodCwAAyg26bvfBRv18yQYdbnXq/NOS9cKMMxTOIScAgJ/gEwldUrinRjNeXK+axhaN6GPTc9ePpdgAAPwKIzc4oTanSxt2H9Lzq3Zp1VcHJEmZCb20eOZ4RVvZhQAA/oVPJnSopc2lj7dX6p9fVGj1Vwfc16kJDbHoqnF9Nefi05Qcyz2+AAD+h3IDD4caW/SHj3fqrc17daip1b08PipcU4al6Y4JA9QvMdrEhAAAHB/lBm4bdtfortc2q7yuWZKUEmvVFWP6aNKwVI3O7M1tEQAAAYFyAxmGoedXfa3ffVgsp8tQTlK0fvPD03X+oGQuvAcACDiUmx6u0dGmX/79M/3riwpJ0hVj+ui304YzURgAELD4BOvBiivqdffSzdpeUa/wUIsWXj5c14zP5N5PAICARrnpARxtTn3y1UFV2JvlaHPpQL1Dy7dVakdVgyQpKcaqP/1srMb1SzA5KQAAp45yE8T21x7Wi/8p0f9t2qeaxpajXg8PtegHp6XokWnDlG7rZUJCAAC6H+UmSDU62nTlH9eown7kzKfUOKtG9Y1XZHiooq2hys1O1AVDUmTrFW5yUgAAuhflJkj9+ZMSVdib1Se+lxb8aJgmDObMJwBAz0C5CUIH6h360+pdkqRfX3K6Jg5NNTkRAAC+w//KB6H/Xv6VmlqcGpUZr0tGpJkdBwAAn2LkJkg0tbRp76HD2lXVoNfWl0mS5k0dwmndAIAeh3ITBLaV23XlH9focKvTvWzi6Sk6KyfRxFQAAJiDchMEXi7YrcOtTkVHhCozIUpZCVF68IdDzY4FAIApKDcBrqmlTe9+Vi5JemHmeEZrAAA9HhOKA9y/v6hQg6NNWQlRys3mCsMAAFBuAtwbhUcmD181ri+ThwEAEOUmoO2pbtTar2tksUjTx/U1Ow4AAH6BchPA/l64V5J03qBkZcRzbygAACQmFAcUp8vQ3wvL9M8vKrS/9rD2VDdKOnJICgAAHEG5CQCGYeiTHQf12L+2aXtFvcdr2UnRupjbKwAA4Ea58VMtbS5tLbfrwy0V+ucX5dpT3SRJiosM020/GKDRmfFKt0UqMyFK4dwQEwAAN8qNn1mxvUq/X75D2/bb1eJ0uZdHhofomvFZuvuiQeodHWFiQgAA/Bvlxo/UNbVq9qub1Nhy5DYKtl7hystJ1KUj03XhkBRFW/nPBQDAifBp6Uf+sna3GlucOi01Rv97wxnKSoji2jUAAHQR5cZPNLc69eJ/dkuS7pwwUP0So80NBABAgGImqp94Y2OZqhtb1Ld3L/1wZLrZcQAACFiUGz/Q5nTpT6u/liTden6Owjj7CQCAk8anqB9Ysma39h46rMToCF01LtPsOAAABDTm3JjIMAwtWrFTv/vwK0nSbT/IUa+IUJNTAQAQ2Ewduenfv78sFovH4/HHHzczktfZm1v1/pcVenVdqea+/pm72My6YIBuOS/H5HQAAAQ+00duFi5cqFtuucX9PDY21sQ03nWosUWXPvuJ9tc1u5dZLNJDPxyqmedkm5gMAIDgYXq5iY2NVVpamtkxvM4wDP2/t7/Q/rpmJcdaNaqvTQnREfrhyAydf1qy2fEAAAgaFsMwDLO+ef/+/dXc3KzW1lZlZWXpuuuu05w5cxQWduzO5XA45HA43M/tdrsyMzNVV1enuLg4X8Q+KW9u2qu5r3+msBCL3rrzHI3oazM7EgAAprHb7bLZbF75/DZ15Oauu+7S2LFjlZCQoDVr1mjevHkqLy/X008/fcz35Ofna8GCBT5MeerKapr00DtbJEn3TBxEsQEAwIu6feTmgQce0BNPPHHcdbZt26YhQ4YctXzx4sW67bbb1NDQIKvV2uF7A23kZu3X1brrtc2qqndoXL/e+tutZ3EdGwBAj+fNkZtuLzcHDhxQdXX1cdfJyclRRMTRd7besmWLhg8fru3bt2vw4MGd+n7e3DinatGKnfqvD4vlMqSBKTFacuN49e0dZXYsAABMF1CHpZKTk5WcfHITZIuKihQSEqKUlJRuTuV7G3fX6KkPiiVJV47to99OG66oCNPnbwMAEPRM+7QtKCjQunXrdMEFFyg2NlYFBQWaM2eOfvrTn6p3795mxeo2u6ubJEl5OYl6+iejzQ0DAEAPYlq5sVqtWrp0qR5++GE5HA5lZ2drzpw5mjt3rlmRulVDc6skKSHm6MNvAADAe0wrN2PHjtXatWvN+vZe1+BokyTFWjkUBQCAL3HajpfUf1NuYig3AAD4FOXGSxqavyk3kZQbAAB8iXLjJY2M3AAAYArKjZe459wwcgMAgE9Rbrykvv2wlDXc5CQAAPQslBsvaR+5Yc4NAAC+Rbnxkgbm3AAAYArKjZe4z5ai3AAA4FOUGy+p57AUAACmoNx4gaPNqZY2lyRGbgAA8DXKjRc0Opzuf1NuAADwLcqNF7TPt4mKCFVoiMXkNAAA9CyUGy+odxy5IzijNgAA+B7lxgvaD0sxmRgAAN+j3HhBwzcjN7GM3AAA4HOUGy+o547gAACYhnLjBVydGAAA81BuvKD9bKloyg0AAD5HufGC9pEb5twAAOB7lBsvYM4NAADmodx4wbdzbsJNTgIAQM9DufGCBkZuAAAwDeXGC5hzAwCAeSg3XsCp4AAAmIdy4wXucsNhKQAAfI5y4wXuOTeM3AAA4HOUGy/gsBQAAOah3HQzl8vgsBQAACai3HSzxpY2978ZuQEAwPcoN92sfdQmPNQiaxibFwAAX+PTt5t9dzKxxWIxOQ0AAD0P5aab1TPfBgAAU1Fuutm3IzfcVwoAADNQbroZt14AAMBclJtuxmngAACYi3LTzbg6MQAA5qLcdLP2kZtoyg0AAKag3HQz95wbDksBAGAKr5WbRx99VGeffbaioqIUHx/f4TqlpaW69NJLFRUVpZSUFP3yl79UW1tbh+sGinoOSwEAYCqvfQK3tLToqquuUl5enl544YWjXnc6nbr00kuVlpamNWvWqLy8XDfccIPCw8P12GOPeSuW13HTTAAAzOW1kZsFCxZozpw5GjFiRIevf/jhh9q6dateeeUVjR49WlOnTtUjjzyiRYsWqaWlxVuxvK6huVUSZ0sBAGAW0+bcFBQUaMSIEUpNTXUvmzx5sux2u7Zs2XLM9zkcDtntdo+HP+E6NwAAmMu0clNRUeFRbCS5n1dUVBzzffn5+bLZbO5HZmamV3N2lXvODSM3AACYokvl5oEHHpDFYjnuY/v27d7KKkmaN2+e6urq3I+ysjKvfr+uYs4NAADm6tIn8L333quZM2ced52cnJxOfa20tDStX7/eY1llZaX7tWOxWq2yWq2d+h5m4FRwAADM1aVP4OTkZCUnJ3fLN87Ly9Ojjz6qqqoqpaSkSJKWLVumuLg4DR06tFu+hy852pxaX1LjvkIxF/EDAMAcXvsELi0tVU1NjUpLS+V0OlVUVCRJGjhwoGJiYjRp0iQNHTpUP/vZz/Tkk0+qoqJCv/nNbzRr1iy/HpnpyOsbyvTwu1vU1OKUJEWEhah3VITJqQAA6Jm8Vm7mz5+vl156yf18zJgxkqQVK1ZowoQJCg0N1Xvvvac77rhDeXl5io6O1owZM7Rw4UJvRfKal9fuVlOLU8mxVl0wOFk/HpepyPBQs2MBANAjWQzDMMwOcSrsdrtsNpvq6uoUFxfn8+/vchka+tD7am51acV9E5SdFO3zDAAABBpvfn5zb6lTVHaoSc2tLlnDQpSVEGV2HAAAejzKzSkqrqiXJA1IjlFoiMXkNAAAgHJzinZUNUiSTkuNMTkJAACQKDen7KvKIyM3g1JjTU4CAAAkys0p+6qyfeSGcgMAgD+g3JwCp8vQrgNHys1gyg0AAH6BcnMK9lQ3qqXNpV7hoerbu5fZcQAAgCg3p6R9vs3AlBiFcKYUAAB+gXJzCtrn2wziTCkAAPwG5eYUtI/cMJkYAAD/Qbk5BTsqmUwMAIC/odycpFanS18f5LAUAAD+hnJzkvZUN6rVaSg6IlR94jlTCgAAfxFmdoBAc6ixRSu/qtK7n5VLkgamxspi4UwpAAD8BeWmCyrtzbr02U91sMHhXnZWToKJiQAAwPdRbrrgfz7eqYMNDmXYIvWj0X107sAknT0g0exYAADgOyg3nVRW06SlG0olSf/1k9HKo9QAAOCXmFDcSX/4eIdanYbOGZhIsQEAwI9Rbjrh6wMN+r9N+yRJ904abHIaAABwPByWOobNpYe0Zle1tlfUq3B3jZwuQxcNSdHYrN5mRwMAAMdBuTmGdz8r1+L/lLifx0WG6b7JjNoAAODvKDfHkDcgUQcbHDo9PU5D0mM1um+8ekdHmB0LAACcAOXmGC4emqqLh6aaHQMAAHQRE4oBAEBQodwAAICgQrkBAABBhXIDAACCCuUGAAAEFcoNAAAIKpQbAAAQVCg3AAAgqFBuAABAUKHcAACAoEK5AQAAQYVyAwAAggrlBgAABBXKDQAACCqUGwAAEFS8Vm4effRRnX322YqKilJ8fHyH61gslqMeS5cu9VYkAADQA4R56wu3tLToqquuUl5enl544YVjrvfiiy9qypQp7ufHKkIAAACd4bVys2DBAknSkiVLjrtefHy80tLSvBUDAAD0MKbPuZk1a5aSkpJ05plnavHixTIM47jrOxwO2e12jwcAAEA7r43cdMbChQt14YUXKioqSh9++KHuvPNONTQ06K677jrme/Lz892jQgAAAN9nMU40VPIdDzzwgJ544onjrrNt2zYNGTLE/XzJkiW65557VFtbe8KvP3/+fL344osqKys75joOh0MOh8P93G63KzMzU3V1dYqLizvxDwEAAExnt9tls9m88vndpZGbe++9VzNnzjzuOjk5OScdJjc3V4888ogcDoesVmuH61it1mO+BgAA0KVyk5ycrOTkZG9lUVFRkXr37k15AQAAJ81rc25KS0tVU1Oj0tJSOZ1OFRUVSZIGDhyomJgYvfvuu6qsrNRZZ52lyMhILVu2TI899pjuu+8+b0UCAAA9gNfKzfz58/XSSy+5n48ZM0aStGLFCk2YMEHh4eFatGiR5syZI8MwNHDgQD399NO65ZZbvBUJAAD0AF2aUOyPvDkhCQAAeIc3P79Nv84NAABAd6LcAACAoEK5AQAAQYVyAwAAggrlBgAABBXKDQAACCqUGwAAEFQoNwAAIKhQbgAAQFCh3AAAgKBCuQEAAEGFcgMAAIIK5QYAAAQVyg0AAAgqlBsAABBUKDcAACCoUG4AAEBQCTM7wKkyDEOSZLfbTU4CAAA6q/1zu/1zvDsFfLmpr6+XJGVmZpqcBAAAdFV9fb1sNlu3fk2L4Y3K5EMul0v79+9XbGysLBZLt35tu92uzMxMlZWVKS4urlu/dqBhW3hie3yLbeGJ7fEttoUntse32rfF1q1bNXjwYIWEdO8smYAfuQkJCVHfvn29+j3i4uJ6/I7Yjm3hie3xLbaFJ7bHt9gWntge3+rTp0+3FxuJCcUAACDIUG4AAEBQodwch9Vq1UMPPSSr1Wp2FNOxLTyxPb7FtvDE9vgW28IT2+Nb3t4WAT+hGAAA4LsYuQEAAEGFcgMAAIIK5QYAAAQVyg0AAAgqlBsAABBUKDfHsGjRIvXv31+RkZHKzc3V+vXrzY7kdfn5+Ro/frxiY2OVkpKiadOmqbi42GOdCRMmyGKxeDxuv/12kxJ718MPP3zUzzpkyBD3683NzZo1a5YSExMVExOj6dOnq7Ky0sTE3tO/f/+jtoXFYtGsWbMkBf9+sXr1al122WXKyMiQxWLR22+/7fG6YRiaP3++0tPT1atXL02cOFE7duzwWKempkbXX3+94uLiFB8fr5tuukkNDQ0+/Cm6x/G2RWtrq+6//36NGDFC0dHRysjI0A033KD9+/d7fI2O9qfHH3/cxz9J9zjRvjFz5syjftYpU6Z4rNMT9g1JHf4NsVgseuqpp9zrdNe+QbnpwN/+9jfNnTtXDz30kDZt2qRRo0Zp8uTJqqqqMjuaV61atUqzZs3S2rVrtWzZMrW2tmrSpElqbGz0WO+WW25ReXm5+/Hkk0+alNj7hg0b5vGzfvrpp+7X5syZo3fffVdvvPGGVq1apf379+vKK680Ma33bNiwwWM7LFu2TJJ01VVXudcJ5v2isbFRo0aN0qJFizp8/cknn9Szzz6r559/XuvWrVN0dLQmT56s5uZm9zrXX3+9tmzZomXLlum9997T6tWrdeutt/rqR+g2x9sWTU1N2rRpkx588EFt2rRJb775poqLi/WjH/3oqHUXLlzosb/84he/8EX8bneifUOSpkyZ4vGzvvbaax6v94R9Q5LHNigvL9fixYtlsVg0ffp0j/W6Zd8wcJQzzzzTmDVrlvu50+k0MjIyjPz8fBNT+V5VVZUhyVi1apV72Q9+8APj7rvvNi+UDz300EPGqFGjOnyttrbWCA8PN9544w33sm3bthmSjIKCAh8lNM/dd99tDBgwwHC5XIZh9Kz9QpLx1ltvuZ+7XC4jLS3NeOqpp9zLamtrDavVarz22muGYRjG1q1bDUnGhg0b3Ov8+9//NiwWi7Fv3z6fZe9u398WHVm/fr0hydizZ497Wb9+/YxnnnnGu+FM0NH2mDFjhnH55Zcf8z09ed+4/PLLjQsvvNBjWXftG4zcfE9LS4sKCws1ceJE97KQkBBNnDhRBQUFJibzvbq6OklSQkKCx/K//vWvSkpK0vDhwzVv3jw1NTWZEc8nduzYoYyMDOXk5Oj6669XaWmpJKmwsFCtra0e+8mQIUOUlZUV9PtJS0uLXnnlFf385z+XxWJxL+9J+8V3lZSUqKKiwmNfsNlsys3Nde8LBQUFio+P1xlnnOFeZ+LEiQoJCdG6det8ntmX6urqZLFYFB8f77H88ccfV2JiosaMGaOnnnpKbW1t5gT0gZUrVyolJUWDBw/WHXfcoerqavdrPXXfqKys1D//+U/ddNNNR73WHftGwN8VvLsdPHhQTqdTqampHstTU1O1fft2k1L5nsvl0j333KNzzjlHw4cPdy+/7rrr1K9fP2VkZOjzzz/X/fffr+LiYr355psmpvWO3NxcLVmyRIMHD1Z5ebkWLFig8847T19++aUqKioUERFx1B/s1NRUVVRUmBPYR95++23V1tZq5syZ7mU9ab/4vvb/3h39zWh/raKiQikpKR6vh4WFKSEhIaj3l+bmZt1///269tprPe6Cfdddd2ns2LFKSEjQmjVrNG/ePJWXl+vpp582Ma13TJkyRVdeeaWys7O1a9cu/frXv9bUqVNVUFCg0NDQHrtvvPTSS4qNjT3qUH537RuUG3Ro1qxZ+vLLLz3mmEjyOA48YsQIpaen66KLLtKuXbs0YMAAX8f0qqlTp7r/PXLkSOXm5qpfv356/fXX1atXLxOTmeuFF17Q1KlTlZGR4V7Wk/YLdE5ra6t+8pOfyDAMPffccx6vzZ071/3vkSNHKiIiQrfddpvy8/OD7r5L11xzjfvfI0aM0MiRIzVgwACtXLlSF110kYnJzLV48WJdf/31ioyM9FjeXfsGh6W+JykpSaGhoUed9VJZWam0tDSTUvnW7Nmz9d5772nFihXq27fvcdfNzc2VJO3cudMX0UwVHx+v0047TTt37lRaWppaWlpUW1vrsU6w7yd79uzRRx99pJtvvvm46/Wk/aL9v/fx/makpaUddUJCW1ubampqgnJ/aS82e/bs0bJlyzxGbTqSm5urtrY27d692zcBTZSTk6OkpCT370ZP2zck6ZNPPlFxcfEJ/45IJ79vUG6+JyIiQuPGjdPy5cvdy1wul5YvX668vDwTk3mfYRiaPXu23nrrLX388cfKzs4+4XuKiookSenp6V5OZ76Ghgbt2rVL6enpGjdunMLDwz32k+LiYpWWlgb1fvLiiy8qJSVFl1566XHX60n7RXZ2ttLS0jz2BbvdrnXr1rn3hby8PNXW1qqwsNC9zscffyyXy+UugsGivdjs2LFDH330kRITE0/4nqKiIoWEhBx1eCYY7d27V9XV1e7fjZ60b7R74YUXNG7cOI0aNeqE6570vnHKU5KD0NKlSw2r1WosWbLE2Lp1q3Hrrbca8fHxRkVFhdnRvOqOO+4wbDabsXLlSqO8vNz9aGpqMgzDMHbu3GksXLjQ2Lhxo1FSUmK88847Rk5OjnH++eebnNw77r33XmPlypVGSUmJ8Z///MeYOHGikZSUZFRVVRmGYRi33367kZWVZXz88cfGxo0bjby8PCMvL8/k1N7jdDqNrKws4/777/dY3hP2i/r6emPz5s3G5s2bDUnG008/bWzevNl9BtDjjz9uxMfHG++8847x+eefG5dffrmRnZ1tHD582P01pkyZYowZM8ZYt26d8emnnxqDBg0yrr32WrN+pJN2vG3R0tJi/OhHPzL69u1rFBUVefwdcTgchmEYxpo1a4xnnnnGKCoqMnbt2mW88sorRnJysnHDDTeY/JOdnONtj/r6euO+++4zCgoKjJKSEuOjjz4yxo4dawwaNMhobm52f42esG+0q6urM6KiooznnnvuqPd3575BuTmGP/zhD0ZWVpYRERFhnHnmmcbatWvNjuR1kjp8vPjii4ZhGEZpaalx/vnnGwkJCYbVajUGDhxo/PKXvzTq6urMDe4lV199tZGenm5EREQYffr0Ma6++mpj586d7tcPHz5s3HnnnUbv3r2NqKgo44orrjDKy8tNTOxdH3zwgSHJKC4u9ljeE/aLFStWdPi7MWPGDMMwjpwO/uCDDxqpqamG1Wo1LrrooqO2U3V1tXHttdcaMTExRlxcnHHjjTca9fX1Jvw0p+Z426KkpOSYf0dWrFhhGIZhFBYWGrm5uYbNZjMiIyON008/3Xjsscc8PuwDyfG2R1NTkzFp0iQjOTnZCA8PN/r162fccsstR/2Pck/YN9r96U9/Mnr16mXU1tYe9f7u3DcshmEYXRvrAQAA8F/MuQEAAEGFcgMAAIIK5QYAAAQVyg0AAAgqlBsAABBUKDcAACCoUG4AAEBQodwAAICgQrkBAABBhXIDAACCCuUGAAAElf8P7Y3kKRa5/nAAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"import matplotlib.pyplot as plt\n",
"\n",
"plt.plot(np.sort(dh.eig()))\n",
"plt.ylim(None, 20)\n",
"np.real(dh.eig()).min()"
]
},
{
"cell_type": "code",
"execution_count": 4,
"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": 5,
"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",
" # if orbital is not set use all\n",
" if \"l\" not in mag_ent.keys():\n",
" mag_ent[\"l\"] = \"all\"\n",
" if isinstance(mag_ent[\"atom\"], int):\n",
" mag_ent[\"tags\"] = [\n",
" f\"[{mag_ent['atom']}]{dh.atoms[mag_ent['atom']].tag}({mag_ent['l']})\"\n",
" ]\n",
" mag_ent[\"xyz\"] = [dh.xyz[mag_ent[\"atom\"]]]\n",
" if isinstance(mag_ent[\"atom\"], list):\n",
" mag_ent[\"tags\"] = []\n",
" mag_ent[\"xyz\"] = []\n",
" # iterate over atoms\n",
" for atom_idx in mag_ent[\"atom\"]:\n",
" mag_ent[\"tags\"].append(\n",
" f\"[{atom_idx}]{dh.atoms[atom_idx].tag}({mag_ent['l']})\"\n",
" )\n",
" mag_ent[\"xyz\"].append(dh.xyz[atom_idx])\n",
"\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",
" pair[\"tags\"] = []\n",
" for mag_ent in [magnetic_entities[pair[\"ai\"]], magnetic_entities[pair[\"aj\"]]]:\n",
" tag = \"\"\n",
" # get atoms of magnetic entity\n",
" atoms_idx = mag_ent[\"atom\"]\n",
" orbitals = mag_ent[\"l\"]\n",
"\n",
" # if magnetic entity contains one atoms\n",
" if isinstance(atoms_idx, int):\n",
" tag += f\"[{atoms_idx}]{dh.atoms[atoms_idx].tag}({orbitals})\"\n",
"\n",
" # if magnetic entity contains more than one atoms\n",
" if isinstance(atoms_idx, list):\n",
" # iterate over atoms\n",
" atom_group = \"{\"\n",
" for atom_idx in atoms_idx:\n",
" atom_group += f\"[{atom_idx}]{dh.atoms[atom_idx].tag}({orbitals})--\"\n",
" # end {} of the atoms in the magnetic entity\n",
" tag += atom_group[:-2] + \"}\"\n",
" pair[\"tags\"].append(tag)\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": 6,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"k loop: 0%| | 0/100 [00:00<?, ?it/s]"
]
}
],
"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\")\n",
"\n",
"times[\"k_set_time\"] = timer()"
]
},
{
"cell_type": "code",
"execution_count": 7,
"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": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Total number of k points: 100\n",
"Total number of directions: 3\n",
"Total number of matrix inversions: 300\n",
"The shape of the Hamiltonian and the Greens function is 84x84=7056\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"k loop: 100%|██████████| 100/100 [00:27<00:00, 3.70it/s]\n"
]
}
],
"source": [
"if rank == root_node:\n",
" print(\"Starting matrix inversions\")\n",
" print(f\"Total number of k points: {kset.shape[0]}\")\n",
" print(f\"Total number of directions: {len(hamiltonians)}\")\n",
" print(f\"Total number of matrix inversions: {kset.shape[0] * len(hamiltonians)}\")\n",
" print(f\"The shape of the Hamiltonian and the Greens function is {NO}x{NO}={NO*NO}\")\n",
" print(\n",
" \"================================================================================================================================================================\"\n",
" )\n",
"\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": 9,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"##################################################################### GROGU OUTPUT #############################################################################\n",
"================================================================================================================================================================\n",
"Input file: \n",
"Not yet specified.\n",
"Number of nodes in the parallel cluster: 1\n",
"================================================================================================================================================================\n",
"Cell [Ang]: \n",
"[[ 3.79100000e+00 0.00000000e+00 0.00000000e+00]\n",
" [-1.89550000e+00 3.28310231e+00 0.00000000e+00]\n",
" [ 1.25954923e-15 2.18160327e-15 2.05700000e+01]]\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",
"Parameters for the contour integral:\n",
"Number of k points: 10\n",
"k point directions: xy\n",
"Ebot: -15\n",
"Eset: 50\n",
"Esetp: 1000\n",
"================================================================================================================================================================\n",
"Atomic information: \n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[atom index]Element(orbitals) x [Ang] y [Ang] z [Ang] Sx Sy Sz Q Lx Ly Lz Jx Jy Jz\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[3]Fe(2) -7.339158738013707e-06 4.149278510690423e-06 11.657585837928032\n",
"\n",
"[4]Fe(2) -7.326987662162937e-06 4.158274523275774e-06 8.912422537596708\n",
"\n",
"[5]Fe(2) 1.8954667088117545 1.0943913231921656 10.285002698393109\n",
"\n",
"================================================================================================================================================================\n",
"Exchange [meV]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Magnetic entity1 Magnetic entity2 [i j k] d [Ang]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[3]Fe(2) [4]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -39.33305904162077\n",
"DMI: [-3.49132045e-01 2.55180589e-04 -2.79328031e-07]\n",
"Symmetric-anisotropy: [-3.25353129e-02 -1.12520210e-01 -2.11441543e-05 2.53477251e-04\n",
" -6.21351953e-07]\n",
"Energies for debugging: \n",
"array([[-3.94455793e-02, -3.49131423e-04, 3.49132666e-04,\n",
" -3.91880035e-02],\n",
" [-3.96260920e-02, -2.53477251e-07, 2.56883927e-07,\n",
" -3.93655944e-02],\n",
" [-3.24859133e-02, 2.08648262e-08, 2.14234823e-08,\n",
" -3.24858874e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.03936559, -0.03248591, -0.03944558])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.03936559435457027 -0.03248588738966994\n",
"\n",
"[4]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -65.55485507687138\n",
"DMI: [-3.57652832e+00 6.16350266e+00 1.76417146e-05]\n",
"Symmetric-anisotropy: [0.08815647 0.19563304 0.03702482 6.25899263 0.04349251]\n",
"Energies for debugging: \n",
"array([[-6.53592220e-02, -3.62002083e-03, 3.53303581e-03,\n",
" -6.58386446e-02],\n",
" [-6.51008759e-02, -6.25899263e-03, 6.06801268e-03,\n",
" -6.54666986e-02],\n",
" [-6.59577470e-02, -3.70071792e-05, -3.70424626e-05,\n",
" -6.59146548e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.0654667 , -0.06595775, -0.06535922])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06546669860888496 -0.06591465478729451\n",
"\n",
"[3]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -65.55361001040538\n",
"DMI: [ 3.55870825e+00 -6.14744383e+00 2.20064459e-05]\n",
"Symmetric-anisotropy: [ 0.08580243 0.19612473 0.03702668 -6.24894432 -0.04979181]\n",
"Energies for debugging: \n",
"array([[-6.53574853e-02, 3.60850006e-03, -3.50891645e-03,\n",
" -6.58355372e-02],\n",
" [-6.51026555e-02, 6.24894432e-03, -6.04594334e-03,\n",
" -6.54678076e-02],\n",
" [-6.59577277e-02, -3.70046728e-05, -3.70486857e-05,\n",
" -6.59146327e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06546781, -0.06595773, -0.06535749])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06546780758377949 -0.06591463268052017\n",
"\n",
"================================================================================================================================================================\n",
"Runtime information: \n",
"Total runtime: 28.020738958 s\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Initial setup: 0.16348762500000014 s\n",
"Hamiltonian conversion and XC field extraction: 0.739 s\n",
"Pair and site datastructure creatrions: 0.016 s\n",
"k set cration and distribution: 0.018 s\n",
"Rotating XC potential: 0.216 s\n",
"Greens function inversion: 26.841 s\n",
"Calculate energies and magnetic components: 0.027 s\n"
]
}
],
"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",
" # calculate magnetic parameters\n",
" for pair in pairs:\n",
" J_iso, J_S, D = calculate_exchange_tensor(pair)\n",
" pair[\"J_iso\"] = J_iso * sisl.unit_convert(\"eV\", \"meV\")\n",
" pair[\"J_S\"] = J_S * sisl.unit_convert(\"eV\", \"meV\")\n",
" pair[\"D\"] = D * sisl.unit_convert(\"eV\", \"meV\")\n",
"\n",
" times[\"end_time\"] = timer()\n",
" print(\n",
" \"##################################################################### GROGU OUTPUT #############################################################################\"\n",
" )\n",
"\n",
" print_parameters(simulation_parameters)\n",
" print_atoms_and_pairs(magnetic_entities, pairs)\n",
" print_runtime_informations(times)\n",
"\n",
" # remove clutter to save magnetic entities and pair information\n",
" for pair in pairs:\n",
" del pair[\"Gij\"]\n",
" del pair[\"Gij_tmp\"]\n",
" del pair[\"Gji\"]\n",
" del pair[\"Gji_tmp\"]\n",
" for mag_ent in magnetic_entities:\n",
" del mag_ent[\"Gii\"]\n",
" del mag_ent[\"Gii_tmp\"]\n",
" del mag_ent[\"Vu1\"]\n",
" del mag_ent[\"Vu2\"]\n",
" results = dict(\n",
" parameters=simulation_parameters,\n",
" magnetic_entities=magnetic_entities,\n",
" pairs=pairs,\n",
" runtime=times,\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'parameters': {'path': 'Not yet specified.',\n",
" 'scf_xcf_orientation': array([0, 0, 1]),\n",
" 'ref_xcf_orientations': [{'o': array([1, 0, 0]),\n",
" 'vw': [array([0, 1, 0]), array([0, 0, 1])]},\n",
" {'o': array([0, 1, 0]), 'vw': [array([1, 0, 0]), array([0, 0, 1])]},\n",
" {'o': array([0, 0, 1]), 'vw': [array([1, 0, 0]), array([0, 1, 0])]}],\n",
" 'kset': 10,\n",
" 'kdirs': 'xy',\n",
" 'ebot': -15,\n",
" 'eset': 50,\n",
" 'esetp': 1000,\n",
" 'parallel_size': 1,\n",
" 'cell': array([[ 3.79100000e+00, 0.00000000e+00, 0.00000000e+00],\n",
" [-1.89550000e+00, 3.28310231e+00, 0.00000000e+00],\n",
" [ 1.25954923e-15, 2.18160327e-15, 2.05700000e+01]])},\n",
" 'magnetic_entities': [{'atom': 3,\n",
" 'l': 2,\n",
" 'orbital_indeces': array([41, 42, 43, 44, 45, 46, 47, 48, 49, 50], dtype=int32),\n",
" 'spin_box_indeces': array([ 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94,\n",
" 95, 96, 97, 98, 99, 100, 101]),\n",
" 'tags': ['[3]Fe(2)'],\n",
" 'xyz': [array([-7.33915874e-06, 4.14927851e-06, 1.16575858e+01])],\n",
" 'energies': [[3.734465086018773, 3.7345045132387904],\n",
" [3.734496010868627, 3.7344992971735755],\n",
" [3.744098645200186, 3.744096314352881]]},\n",
" {'atom': 4,\n",
" 'l': 2,\n",
" 'orbital_indeces': array([56, 57, 58, 59, 60, 61, 62, 63, 64, 65], dtype=int32),\n",
" 'spin_box_indeces': array([112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124,\n",
" 125, 126, 127, 128, 129, 130, 131]),\n",
" 'tags': ['[4]Fe(2)'],\n",
" 'xyz': [array([-7.32698766e-06, 4.15827452e-06, 8.91242254e+00])],\n",
" 'energies': [[3.734489896651657, 3.7345163491510975],\n",
" [3.7346033722527983, 3.73451361626319],\n",
" [3.744120742400213, 3.744123051274557]]},\n",
" {'atom': 5,\n",
" 'l': 2,\n",
" 'orbital_indeces': array([71, 72, 73, 74, 75, 76, 77, 78, 79, 80], dtype=int32),\n",
" 'spin_box_indeces': array([142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154,\n",
" 155, 156, 157, 158, 159, 160, 161]),\n",
" 'tags': ['[5]Fe(2)'],\n",
" 'xyz': [array([ 1.89546671, 1.09439132, 10.2850027 ])],\n",
" 'energies': [[2.1707567909611734, 2.170785438112471],\n",
" [2.1708196859983695, 2.1707887369208603],\n",
" [2.171678783995076, 2.1716787861431013]]}],\n",
" 'pairs': [{'ai': 0,\n",
" 'aj': 1,\n",
" 'Ruc': array([0, 0, 0]),\n",
" 'tags': ['[3]Fe(2)', '[4]Fe(2)'],\n",
" 'energies': [[-0.03944557925189693,\n",
" -0.0003491314233520734,\n",
" 0.00034913266605597995,\n",
" -0.03918800351839513],\n",
" [-0.03962609204615879,\n",
" -2.5347725119651397e-07,\n",
" 2.5688392741602586e-07,\n",
" -0.03936559435457027],\n",
" [-0.032485913309031784,\n",
" 2.086482622369864e-08,\n",
" 2.1423482286665388e-08,\n",
" -0.03248588738966994]],\n",
" 'J_iso': -39.33305904162077,\n",
" 'J_S': array([-3.25353129e-02, -1.12520210e-01, -2.11441543e-05, 2.53477251e-04,\n",
" -6.21351953e-07]),\n",
" 'D': array([-3.49132045e-01, 2.55180589e-04, -2.79328031e-07])},\n",
" {'ai': 1,\n",
" 'aj': 2,\n",
" 'Ruc': array([0, 0, 0]),\n",
" 'tags': ['[4]Fe(2)', '[5]Fe(2)'],\n",
" 'energies': [[-0.06535922203590024,\n",
" -0.003620020827251614,\n",
" 0.0035330358091578322,\n",
" -0.06583864458582893],\n",
" [-0.065100875948478,\n",
" -0.006258992627124543,\n",
" 0.0060680126832299295,\n",
" -0.06546669860888496],\n",
" [-0.06595774697036409,\n",
" -3.700717915377913e-05,\n",
" -3.7042462582999624e-05,\n",
" -0.06591465478729451]],\n",
" 'J_iso': -65.55485507687138,\n",
" 'J_S': array([0.08815647, 0.19563304, 0.03702482, 6.25899263, 0.04349251]),\n",
" 'D': array([-3.57652832e+00, 6.16350266e+00, 1.76417146e-05])},\n",
" {'ai': 0,\n",
" 'aj': 2,\n",
" 'Ruc': array([0, 0, 0]),\n",
" 'tags': ['[3]Fe(2)', '[5]Fe(2)'],\n",
" 'energies': [[-0.0653574852836104,\n",
" 0.003608500059552775,\n",
" -0.0035089164454806124,\n",
" -0.06583553716382626],\n",
" [-0.06510265546673727,\n",
" 0.006248944319269189,\n",
" -0.006045943341346723,\n",
" -0.06546780758377949],\n",
" [-0.0659577276970226,\n",
" -3.700467282282758e-05,\n",
" -3.704868571469545e-05,\n",
" -0.06591463268052017]],\n",
" 'J_iso': -65.55361001040538,\n",
" 'J_S': array([ 0.08580243, 0.19612473, 0.03702668, -6.24894432, -0.04979181]),\n",
" 'D': array([ 3.55870825e+00, -6.14744383e+00, 2.20064459e-05])}],\n",
" 'runtime': {'start_time': 3.2424505,\n",
" 'setup_time': 3.405938125,\n",
" 'H_and_XCF_time': 4.145158916,\n",
" 'site_and_pair_dictionaries_time': 4.161620875,\n",
" 'k_set_time': 4.179151,\n",
" 'reference_rotations_time': 4.395303708,\n",
" 'green_function_inversion_time': 31.236585666,\n",
" 'end_time': 31.263189458}}"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"results"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"ename": "SyntaxError",
"evalue": "invalid syntax (4136980242.py, line 2)",
"output_type": "error",
"traceback": [
"\u001b[0;36m Cell \u001b[0;32mIn[11], line 2\u001b[0;36m\u001b[0m\n\u001b[0;31m ================================================================================================================================================================\u001b[0m\n\u001b[0m ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m invalid syntax\n"
]
}
],
"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",
"[[ 3.79100000e+00 0.00000000e+00 0.00000000e+00]\n",
" [-1.89550000e+00 3.28310231e+00 0.00000000e+00]\n",
" [ 1.25954923e-15 2.18160327e-15 2.05700000e+01]]\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: 20\n",
"k point directions: xy\n",
"================================================================================================================================================================\n",
"Parameters for the contour integral:\n",
"Ebot: -15\n",
"Eset: 50\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",
"[3]Fe(2) -7.339158738013707e-06 4.149278510690423e-06 11.657585837928032\n",
"\n",
"[4]Fe(2) -7.326987662162937e-06 4.158274523275774e-06 8.912422537596708\n",
"\n",
"[5]Fe(2) 1.8954667088117545 1.0943913231921656 10.285002698393109\n",
"\n",
"================================================================================================================================================================\n",
"Exchange [meV]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Magnetic entity1 Magnetic entity2 [i j k] d [Ang]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[3]Fe(2) [4]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -63.47396322818869\n",
"DMI: [-9.32946318e-01 -7.94980656e-04 -1.68626155e-06]\n",
"Symmetric-anisotropy: [-3.28657093e+00 1.26919731e+00 -5.25251097e-04 -7.51662131e-04\n",
" -1.69771318e-05]\n",
"Energies for debugging: \n",
"array([[-6.22047659e-02, -9.32929341e-04, 9.32963295e-04,\n",
" -6.14565896e-02],\n",
" [-6.75292273e-02, 7.51662131e-07, -8.38299181e-07,\n",
" -6.67605342e-02],\n",
" [-5.46814050e-02, 5.23564835e-07, 5.26937358e-07,\n",
" -5.46808408e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06676053, -0.0546814 , -0.06220477])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06676053415896452 -0.05468084079537627\n",
"\n",
"[4]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.97346468960899\n",
"DMI: [-3.79957858e+00 6.15341580e+00 1.98171298e-03]\n",
"Symmetric-anisotropy: [0.10023044 0.10625949 0.07113812 6.2430436 0.03640581]\n",
"Energies for debugging: \n",
"array([[-6.08672052e-02, -3.83598439e-03, 3.76317277e-03,\n",
" -6.11799546e-02],\n",
" [-6.06363005e-02, -6.24304360e-03, 6.06378800e-03,\n",
" -6.08732343e-02],\n",
" [-5.98939156e-02, -6.91564033e-05, -7.31198293e-05,\n",
" -5.98115156e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06087323, -0.05989392, -0.06086721])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06087323425105703 -0.05981151560602731\n",
"\n",
"[3]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.97783221534475\n",
"DMI: [ 3.78495282e+00 -6.14017205e+00 1.98157781e-03]\n",
"Symmetric-anisotropy: [ 0.09115149 0.10676698 0.07113987 -6.23329445 -0.04250828]\n",
"Energies for debugging: \n",
"array([[-6.08710652e-02, 3.82746110e-03, -3.74244454e-03,\n",
" -6.11757507e-02],\n",
" [-6.06489192e-02, 6.23329445e-03, -6.04704964e-03,\n",
" -6.08866807e-02],\n",
" [-5.98938645e-02, -6.91582962e-05, -7.31214518e-05,\n",
" -5.98114592e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06088668, -0.05989386, -0.06087107])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06088668072635317 -0.05981145916679715\n",
"\n",
"================================================================================================================================================================\n",
"Runtime information: \n",
"Total runtime: 106.24268275 s\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Initial setup: 0.12151341699999985 s\n",
"Hamiltonian conversion and XC field extraction: 0.664 s\n",
"Pair and site datastructure creatrions: 0.008 s\n",
"k set cration and distribution: 0.013 s\n",
"Rotating XC potential: 0.209 s\n",
"Greens function inversion: 105.202 s\n",
"Calculate energies and magnetic components: 0.025 s"
]
},
{
"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",
"[[ 3.79100000e+00 0.00000000e+00 0.00000000e+00]\n",
" [-1.89550000e+00 3.28310231e+00 0.00000000e+00]\n",
" [ 1.25954923e-15 2.18160327e-15 2.05700000e+01]]\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: 20\n",
"k point directions: xy\n",
"================================================================================================================================================================\n",
"Parameters for the contour integral:\n",
"Ebot: -15\n",
"Eset: 50\n",
"Esetp: 1000\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",
"[3]Fe(2) -7.339158738013707e-06 4.149278510690423e-06 11.657585837928032\n",
"\n",
"[4]Fe(2) -7.326987662162937e-06 4.158274523275774e-06 8.912422537596708\n",
"\n",
"[5]Fe(2) 1.8954667088117545 1.0943913231921656 10.285002698393109\n",
"\n",
"================================================================================================================================================================\n",
"Exchange [meV]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Magnetic entity1 Magnetic entity2 [i j k] d [Ang]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[3]Fe(2) [4]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -63.613204102235464\n",
"DMI: [-9.33066197e-01 -9.12173634e-04 -2.66552375e-06]\n",
"Symmetric-anisotropy: [-3.41993761e+00 1.35139766e+00 -4.15748909e-03 -8.08118347e-04\n",
" 5.75970072e-05]\n",
"Energies for debugging: \n",
"array([[-6.22618064e-02, -9.33123794e-04, 9.33008600e-04,\n",
" -6.15446642e-02],\n",
" [-6.76703372e-02, 8.08118347e-07, -1.01622892e-06,\n",
" -6.70331417e-02],\n",
" [-5.07756491e-02, 4.15482357e-06, 4.16015462e-06,\n",
" -5.07708822e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06703314, -0.05077565, -0.06226181])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06703314171283804 -0.05077088224506824\n",
"\n",
"[4]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.981389719594276\n",
"DMI: [-3.79876136e+00 6.14996306e+00 4.52338839e-03]\n",
"Symmetric-anisotropy: [0.07285133 0.12581883 0.07081179 6.23894685 0.0361525 ]\n",
"Energies for debugging: \n",
"array([[-6.08555709e-02, -3.83491386e-03, 3.76260887e-03,\n",
" -6.11800599e-02],\n",
" [-6.06231936e-02, -6.23894685e-03, 6.06097927e-03,\n",
" -6.09085384e-02],\n",
" [-5.99055239e-02, -6.62883987e-05, -7.53351755e-05,\n",
" -5.98235127e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06090854, -0.05990552, -0.06085557])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06090853838909634 -0.05982351269907684\n",
"\n",
"[3]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.99043921106905\n",
"DMI: [ 3.78566701e+00 -6.13454413e+00 4.52301283e-03]\n",
"Symmetric-anisotropy: [ 0.06072775 0.1281154 0.07081275 -6.22754151 -0.04243893]\n",
"Energies for debugging: \n",
"array([[-6.08623238e-02, 3.82810594e-03, -3.74322808e-03,\n",
" -6.11792824e-02],\n",
" [-6.06338592e-02, 6.22754151e-03, -6.04154675e-03,\n",
" -6.09297115e-02],\n",
" [-5.99054929e-02, -6.62897359e-05, -7.53357616e-05,\n",
" -5.98234762e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06092971, -0.05990549, -0.06086232])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06092971146100409 -0.05982347618612315\n",
"\n",
"================================================================================================================================================================\n",
"Runtime information: \n",
"Total runtime: 107.15435087500009 s\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Initial setup: 0.11759866600004898 s\n",
"Hamiltonian conversion and XC field extraction: 0.737 s\n",
"Pair and site datastructure creatrions: 0.012 s\n",
"k set cration and distribution: 0.005 s\n",
"Rotating XC potential: 0.241 s\n",
"Greens function inversion: 106.014 s\n",
"Calculate energies and magnetic components: 0.028 s"
]
},
{
"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",
"[[ 3.79100000e+00 0.00000000e+00 0.00000000e+00]\n",
" [-1.89550000e+00 3.28310231e+00 0.00000000e+00]\n",
" [ 1.25954923e-15 2.18160327e-15 2.05700000e+01]]\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: 20\n",
"k point directions: xy\n",
"================================================================================================================================================================\n",
"Parameters for the contour integral:\n",
"Ebot: -15\n",
"Eset: 500\n",
"Esetp: 1000\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",
"[3]Fe(2) -7.339158738013707e-06 4.149278510690423e-06 11.657585837928032\n",
"\n",
"[4]Fe(2) -7.326987662162937e-06 4.158274523275774e-06 8.912422537596708\n",
"\n",
"[5]Fe(2) 1.8954667088117545 1.0943913231921656 10.285002698393109\n",
"\n",
"================================================================================================================================================================\n",
"Exchange [meV]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Magnetic entity1 Magnetic entity2 [i j k] d [Ang]\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"[3]Fe(2) [4]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -63.440713197129476\n",
"DMI: [-9.32927095e-01 -6.78818440e-04 -1.30891264e-06]\n",
"Symmetric-anisotropy: [-3.23494810e+00 1.24125542e+00 -3.51321370e-05 -6.74519122e-04\n",
" -2.75924032e-05]\n",
"Energies for debugging: \n",
"array([[-6.21994578e-02, -9.32899503e-04, 9.32954688e-04,\n",
" -6.14470205e-02],\n",
" [-6.74870164e-02, 6.74519122e-07, -6.83117758e-07,\n",
" -6.66756613e-02],\n",
" [-5.57434879e-02, 3.38232244e-08, 3.64410497e-08,\n",
" -5.57434496e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06667566, -0.05574349, -0.06219946])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06667566129701633 -0.05574344955698773\n",
"\n",
"[4]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.970992244504004\n",
"DMI: [-3.79968813e+00 6.15433793e+00 2.99893936e-04]\n",
"Symmetric-anisotropy: [0.10892357 0.1011052 0.07120461 6.24418495 0.03644165]\n",
"Energies for debugging: \n",
"array([[-6.08698870e-02, -3.83612979e-03, 3.76324648e-03,\n",
" -6.11810210e-02],\n",
" [-6.06394618e-02, -6.24418495e-03, 6.06449091e-03,\n",
" -6.08620687e-02],\n",
" [-5.98919690e-02, -7.09047183e-05, -7.15045062e-05,\n",
" -5.98094909e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06086207, -0.05989197, -0.06086989])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.06086206867258083 -0.059809490859558875\n",
"\n",
"[3]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n",
"Isotropic: -60.974228876874314\n",
"DMI: [ 3.78485228e+00 -6.14192269e+00 2.99131471e-04]\n",
"Symmetric-anisotropy: [ 0.10123774 0.10087947 0.07120733 -6.23504361 -0.04251795]\n",
"Energies for debugging: \n",
"array([[-6.08733494e-02, 3.82737022e-03, -3.74233433e-03,\n",
" -6.11763461e-02],\n",
" [-6.06540016e-02, 6.23504361e-03, -6.04880177e-03,\n",
" -6.08729911e-02],\n",
" [-5.98919105e-02, -7.09081946e-05, -7.15064575e-05,\n",
" -5.98094279e-02]])\n",
"J_ii for debugging: (check if this is the same as in calculate_exchange_tensor)\n",
"array([-0.06087299, -0.05989191, -0.06087335])\n",
"Test J_xx = E(y,z) = E(z,y)\n",
"-0.060872991134857855 -0.05980942785099304\n",
"\n",
"================================================================================================================================================================\n",
"Runtime information: \n",
"Total runtime: 995.791283916 s\n",
"----------------------------------------------------------------------------------------------------------------------------------------------------------------\n",
"Initial setup: 0.10009837499999819 s\n",
"Hamiltonian conversion and XC field extraction: 0.522 s\n",
"Pair and site datastructure creatrions: 0.025 s\n",
"k set cration and distribution: 0.023 s\n",
"Rotating XC potential: 0.239 s\n",
"Greens function inversion: 994.703 s\n",
"Calculate energies and magnetic components: 0.178 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 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\")"
]
}
],
"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
}