From cd6e293d23adf45aba4a50319bcce0eef028d86e Mon Sep 17 00:00:00 2001 From: Daniel Pozsar Date: Tue, 29 Oct 2024 14:17:57 +0100 Subject: [PATCH] a bunch of things and proper documentation and build --- .gitignore | 3 + README.md | 65 ++- .../source/_build/doctrees/environment.pickle | Bin 27655 -> 36421 bytes docs/source/_build/doctrees/grogu.doctree | Bin 23697 -> 33107 bytes docs/source/_build/doctrees/index.doctree | Bin 6221 -> 6050 bytes docs/source/_build/doctrees/modules.doctree | Bin 2904 -> 2745 bytes docs/source/_build/html/.buildinfo | 4 +- .../_build/html/_modules/grogu/useful.html | 240 +++++--- docs/source/_build/html/_modules/index.html | 2 +- docs/source/_build/html/_static/basic.css | 15 +- docs/source/_build/html/_static/doctools.js | 7 + .../_build/html/_static/language_data.js | 7 + .../source/_build/html/_static/searchtools.js | 38 +- docs/source/_build/html/genindex.html | 69 ++- docs/source/_build/html/grogu.html | 132 +++-- docs/source/_build/html/index.html | 3 +- docs/source/_build/html/modules.html | 15 +- docs/source/_build/html/objects.inv | 6 +- docs/source/_build/html/py-modindex.html | 2 +- docs/source/_build/html/search.html | 2 +- docs/source/_build/html/searchindex.js | 2 +- local_operator.ipynb | 2 +- pyproject.toml | 30 +- requirements-dev.txt | 1 - requirements.txt | 1 - src/{grogu => grogu_magn}/__init__.py | 0 src/{grogu => grogu_magn}/example.py | 0 src/grogu_magn/grogu.py | 526 ++++++++++++++++++ src/{grogu => grogu_magn}/jij.py | 2 +- src/{grogu => grogu_magn}/useful.py | 0 test.ipynb | 160 ++++-- test.py | 515 ----------------- 32 files changed, 1107 insertions(+), 742 deletions(-) rename src/{grogu => grogu_magn}/__init__.py (100%) rename src/{grogu => grogu_magn}/example.py (100%) create mode 100644 src/grogu_magn/grogu.py rename src/{grogu => grogu_magn}/jij.py (99%) rename src/{grogu => grogu_magn}/useful.py (100%) delete mode 100644 test.py diff --git a/.gitignore b/.gitignore index deae2a3..b3eaab6 100644 --- a/.gitignore +++ b/.gitignore @@ -24,3 +24,6 @@ tmp* # Mac stuff *.DS_Store* + +# PYPI test tokens +.pypirc diff --git a/README.md b/README.md index fb52fed..f9e167c 100644 --- a/README.md +++ b/README.md @@ -1,34 +1,71 @@ # Relativistic magnetic interactions from non-orthogonal basis sets +More on the theoretical background can be seen on [arXiv](https://arxiv.org/abs/2309.02558). + # TODO -- Definition of magnetic entities: - * Through simple sequence o forbitals in the unit cell - * Through atom specification - * Through atom and orbital specification -- Separation of TR and TRB components of the Hamiltonian, Identification of the exchange field. -- Definition of commutator expressions, old projection matrix elements -- Efficient calculation of Green's functions -- Calculation of energy and momentum resolved derivatives -- Parallel BZ and serial energy integral +- Parallel or serial energy integral to reduce memory overhead +- Check the symmetrization of the Hamiltonian and overlap matrix to make them hermitian +- Check if exchange field has scalar part +- Add more tests +- Run tests on different magnetic materials and compare it to Grogu Matlab. +- Ehm MAKE IT WORK SOMEHOW :'( # Building wheel See detailed documentation on [PYPI](https://packaging.python.org/en/latest/tutorials/packaging-projects/). +- First you need some API Tokens for Test PYPI ,to be able to upload. You can read about it [here](https://test.pypi.org/help/#apitoken). I own the current project, so you have to contact me. + + Use the following commands for a quick setup: - Build wheel + ``` python -m build +``` + +- Push to PYPI test repository +``` +python -m twine upload --repository testpypi dist/* ``` -Build wheel: +# Usage -Push to pypi(testpypi for beginners): python3 -m twine upload --repository testpypi dist/* +## For end users +Download and install from [PYPI](https://test.pypi.org/project/grogu-magn/) test repository. -Ehhez kellenek tokenek: -You will be prompted for a username and password. For the username, use __token__. For the password, use the token value, including the pypi- prefix. +``` +python3 -m pip install --index-url https://test.pypi.org/simple/ grogu_magn +``` +## For developers -Végfelhasználóknak (egyelőre testpypi): python3 -m pip install --index-url https://test.pypi.org/simple/ example-package-YOUR-USERNAME-HERE +- Clone repository from Gitea + +``` +git clone https://gitea.vo.elte.hu/et209d/grogu.git +``` + +- Create .venv (for example with VsCode) + * Use python 3.9.6 + * install dependencies from: + * requirements.txt + * requirements-dev.txt + * /docs/requirements.txt + +- Install and run pre-commit + +``` +pre-commit install +pre-commit run --all-files +``` + +To build the documentation navigate to the `docs/source` folder and run `make clean` and `make html`. After this the html page can be found in `docs/source/_build/html`. + +``` +cd docs/source +make clean +make html +``` diff --git a/docs/source/_build/doctrees/environment.pickle b/docs/source/_build/doctrees/environment.pickle index 74e684844b41ad926b8d447fb0873dd23558cd98..2b1ac65c421cc03d651ff95a0f02d0b49d7253e0 100644 GIT binary patch literal 36421 zcmeHwdwg8Sbsk?J9wbPBZ%VRs1?sUNfdv6leo%rehy+Dw34jJbO9}+9@9tjU-izJ6 ztG#;(fHED+4~bB&EyY^Paq2vqHg(-cS~qUoBz2NHb?UT9uR6|3wN z3r4Y0^}8PStbA@X0(z<1`%Lv(HQ{g9OQq_B-;r~&m9k;vta8;qAaQsVE1yi4^op6c zGF87VJzXhZLazD}8g3n^`Yo4TqwH4w-8sEz8TpcP)z#6y{hYH@%sYAxMJ)?uXThub zJJ12LJ6*>sXALCM8KMbTWzTN4Qk|j->@@e7d(AFUnstiv)Ky0}7_?kT>R!b`_f>zt;u?TgEMqC_7HAjVuJ#@`&KtPEcq%VCfwfb@yuC(I(3de)Zx?pLULU zu3eH$@0N|zS*Ns8wiZmF*)?=>=yabZk+p{<=Ahqt1)>92s{R&`TsHEC?qWnOCzB_W zgZKq!^p&c4+HdDndcIHriATVFuWY0h3)1se-YDt?v^?~5*>Nh@l25C$Yfppq*Q&m@ zU{oqdruBT@SxV1)SFf(5D`mat=Jg88W!1kel22nK3uVJ~)7iZ4x`ta_^*ai>TQSOM zmR5^Z%o)q11guPH*PEZWmZ8r6j%3*-eMjp`*~?ZyYL2?Dx_+0joXvYVgY+4ggf^Jn ze*28+6*6uz2Mupeo|}e}J7oi6?y`z?4g77#k3+1xbH=>x98V!0mz=UCcM>q3E;i?Kr1)GtF)7rk7WsMS_jnN-IznF$*JR0t4M$(2FHZaX0Oi^FZHLv??VlySPHU zl(|)c0J6x zMycv=O`_jQ-c7^6So14XPPr@;vs`gWR3281G?EU|V2-6GErc&}LJ2Ef-cH5Bf|e zpQ9F-KSs&R&KJrkTTwXgEj>%0HZwFZUJQma(lQUCIK5>9d}O_(BJ*8#=?^%9`; zl9eOv-XiJcsf~$5N^BKL=U_^hM5sbqSaVoPMv^n!Y}sOyL(C9BEc2Mr87$ZE3|{uT zx#m*YBHbqA3Va<{C#|cfr{^g{m3@cwH?Uiz{wOCqs z?G@-?#8$1W8Ip zz5>`MnBsVql1D4Fa^70ldBcon{8*8&+G6>HjU`SBX$bpH22{Ew{|c;-)~kbCHKt0KsHxQ+HTIsTUW}-qEQ(rmI?#1 z?n9>r4EMePaI?Qq$oFT75>&2=1%{ACM!H9JAT0x#>upaa;9Ni13pb3{Bx;hgpoeft zMfoh1&!IeM4)M`)Bv&!CI0c_Ufs_j>Wk!0SK;*6SY70G3Fpo^}s=rU=X`>)FEZ6;8 z0%$n#DtZ}Ew!Mc{F;8LT*OC)|^!(+E7qy9z@lox8(eaV7iw}*QJ2y2tJ^di2a@b02 zF5p(_1?bJik%>nxk351V7cV}Pw~|hI0o8Yr6Sl0(*-or(9@{+umiDPjvf=(1|S3tfX3^P`#MNLP+N5 zNQhPfBZJ)bwmY%>Xlqp=SC2j4UjMeRP8m7}Gg81}ph$Y;q4KYYz(b^-*ynjgvHh^F z8aa@2EP}@IHqI@Sf{_;swASpZEAOgxv7ea4*d3tZsnm29(q5b3bNRiyWuvB*iSes zisr0IE4JTye7OIeJty9dtt*!njU{jeE3?1DD`HpU?`JUW7=vYmR~Cg!y5_N!rDGYc z4m`q2g_gIy*|e2QV~fk~5>P|^_L1)a6JG=F%BF?qv1qJ7nGmg^iPZ|v0J{UCmQ&6= zsKaD~SJg_^GfzA3cJlnJBCV*{8ewt5W&&<)r}T7a*Ivf*zPn=)3#UE?GEay0>r68#9$_EsA+y6R78Bl z1mKMhxiZYgxO2)5Bg4teDc289>r+^9z;urMaBK!1mW0XeNlOu~`xUJ8QE+tW>L%}l zLt*65Aef1n)xvo+qB1N`>4L$Vc1&Y?2Gnm4Tl#mYhCw+dVzBv3 z3zw>j2sHu`VJ|#MFk|58@^X8m9*c}q_$@pa;9aOHGP;EQKHP;etmQ4L#44JIR#Ygy zD27O8_oynRXq8dEf)KBW`;jK>RdHfDWSUI@2q95cABqpG`dR>=cQ^zpn1mTGq{b#b zs6#7)jKNSpM^)rGz2GHL# znn`oO{^-^%TaY}7$o>wpQ3M|{2koB_@Kfe#RpmbO^(y&>5c@B=~B*`WE8=A5eewg8wm z->%x|0hO5`l?`ainTBe8V!}KPR?V9W{A&iatyKot<|6;*%>w@xgWjCfR)wTwKBbaZ z0|ZY8HQy1yt_38#)BHiz|^#b z!7m0OQuBit=ntD8uToU>Ze+kL?qWLjZ@4p6<@K1tqeaZZ9 zs`$(1S5)$==6_eoZK?d`HzWzkA*Om{=12N09oTs|xnLK&Stg z`DaRfzi<9;RsQGZU&y4tTSjq;*at$B+DbV*X2U)ep-0T_y`oepB5y6t2mN+F4UnX2 z!87hxe6(%;Wu>Y!EB{KE(Z4qThO+h+4212c>W-!s!T~u%2DrjjRNyA|Hl(gIY&%l+ z4*cO<3lgI4PJ*{8cpK%~Rjz|_yJXJ%KjvR@+dTx`tLk)7u3P5(J;C^eqkfIqV(+6$ z`&E+zlsl+$hbVVg<&IG9sLI_!xm#84HY8G8nTFR<@^((%PRSNd-hreTt%j6+C;l+} z7!tjlx(lfwcEgan0TFHPQ7ONjVlU#E0zw0jpE}L4D`IUBH82PC)K)1Z`(88`EKR5s zPEugkfc+;N6YL%21B}23aIA!*IJNpB`@@nkV3)UT+2Y?ZSu`kWWoVJDg67O?^7I2@ z5aX7cHS$a0HtRmfG3fjkfqoS2H_0E86f){4tv0-V72TQM_3f5xky9% zZP82`nN?No?U52c5-x!wUlP$V_{D-aoaqR<=g=0MnxJfNSSI`(ze7xx3y1?$6^2{G zS}Fj6hzva{Mieg?W+eU+)Dnjo<`gDd2iV|45ZaC{Vv3LI3f}8;ON)~ z10fb1{7?~&poZVZW}7yRJbUJWP%S9-LL}qwl`)K;>o ztzvQJ2p=^RN3dQ%pj(|HYZZsH5S)^$&5M2)M+p)1qa7?;uWM)^ofkS`W`Bbgqjq-V z`N*05H3ICCs08Ti0C~VAaS%EVR00>bx4(r7+|d3urKrCB9ctud%v02(Jh_k*4=mhd zwwl`RACbHWfF_P3B0h!ZLtFjBwS|_5GEkX1`4}|B$eD9H;yUN;jgqh{pI;> z2=>y1VOjAkdA_yd){o!w(Ytr@OnJWb^T$PwX3O)fxk1i(bR?6Wa^tK99j5em;hAt= zMc@2gGH#29;1wKR7H7kt5~#SHCeBOM_u=d@cSn}tXfio~|8Np{fEPYID=_B0JY-n) zZ>cSl(M>KqR9k3UHW%DlTaZQw0KwQam!1=)N!Yt912}3uAe7})@=XISp5yf6xxpeG z>=|&gR%vBG#$N~02(+6xiG=49gyu!pZ}j>@3#+*X0@X_OiYtj5XQ)wycSR;doq0%JTrDpP^$pkXt8|vnxGQwCth_s;G9zhbD`iJoYVn^}gck6DhCIBPov>!#W{W$)BMmvRHS_sDR zi>8y}=PlxATKrt*AJ2Xazls80B-W>3I}}&)M)YZ*{SEZ#`5k!&CqyaG6pEF(Ef_*$ z5~~;WIUs$7B0;lT5I}R|uQlq-nSu6_o zFO>C?xkh6_OnYE4*Q&;%Y^ZL6rPrHgDNRoCFcG$LCkXHm)faw5hj{YInr5Ds@_v?> z$}007V(Mm8=9vJ+HAZ1C;ir-2ChT{tQDVP)t%B#)Dqz-iOnVvSS4GwTXVI!TRoo!o z{9VVeo~bVDc=$8u_Xc_RQ=3KcG#Uy?{y3v})tP;O&0zn0T$oQ{r7R81c zpV+LO*)KRR9JEHv@7nQNY5svY^92d@{O<$-H)we}7{&h^>&&IzeY~&FLU$&iiqz`3kq?|#c=2nwBDW{ue#zv)luwKd-LA{`# zY@UJ*3i^qrG2e7ScjFa-CA@BeEhWy+<(6v%{=$t5+;5sG8x{EL;sUQ3Nxk4#o2O)h zf`7JY%r{-|x8i(Y)}x1&2)yEyG9rmMZdr7la*Yc7(v7RYr6R|;Yf3a!I zH(lIq{MJIP3wNmHh8Oq!K+&AbyHQaes1Y^2cG19z>uR2W4T`z5Y0NiW%q^z7SYyo3 z-nfh(YMK!nm2r1m#-P4_9uGI)`vwJkrfJMKUBGw)5^qP*s{we=wnoyW8<%vUX^w1E z(znDV9UVr!$nECI*r3Q4n#O$7MUD#?Ao4{#yQ-1*kKMSuf3#_KY*gNd8_1i+P%rEc zH&4R`gVH!(P#6|HLY12Rg2Epqq2I1c<*M#7Cv2hX0s@2 zdATq>@TAeTy=*#BINaL+i=KD7ni5`?PujIKZ6@Q7V-|gioVmdayTT3RHwuh z0NTXmA~?&5d)qq1X)mXkr(0G!a6m64inFR0`y0eDpU%q~)zlud%YW0rWgHH52cmbv z4g@#C$`h4z#B)HL`5TCwqUQ5XzsEAz#1Xl09tttu&^-JkQ}Jp(YvqGc`^JRpuIa^` z#+Mqou_x$lt<8FJjg4^KeDU`mle>m~bloqzmdjsw*sKNpH|GXQ6yY4&wpkRR5fsed zg#4rt?qSmPoL66Pg@Oh7&hJg_rUCy;!g4TuN}) zFdaOmH&{~p`2cTf8=u}yZT0*;R24^A#9Bb7eFX?N>qS6aasfmfg0Dyxr-I4k;3zY% zNsJH~mLpz*;T4PCfGOjuCA#;UNT8ZHZ_hNm$W#K_yP)drI4YUtLomFNv2iN&S_Ccokr({7%ybF7VaA;E{tSw z4=YquL-_}C{?-t-#Ij%;do0G6!V`bx+I>DNV+XyXp%81s%wIykF#66u^ z;j4j&&j(^qg@y!x+d}ZBg=&M-`*@G4Ra_rUoHtLUumZP>r}q?2rSZ4QY*5~oRj&e} zn3|x4gs_-9RgPXKA&BR1$30kt%ziI4(7&bT#U2y}T(qA>jW`GG=K-3qe~A9#O(}Y} zak4x~H?r?hcSK*JyIJURl-&w4QpQ_46ZZSL=?C~P&F2aG1*BK~_K*oUxZlZqIUBw$ z0PmL$&`XY(56XesWM*g>r1EpsM8e-u(&?W2gJQMLWFx#4ktK*)i=|lTX=s#& zd!rA{)@9ex*pQ;lVhHbwJXu?H#!>Ng2){Gl zL!&5kbPV@IqLu4Xek6uB(kezTJbsSotwZ3x7$OltlN7P(>5Q44D(>fHQzXvci_Q#Q zQL*~cQVO5Q^S9yhU9x~L!x~^k2{x5#@z}v{Z$@f#;#%}Ht}2yRR4MDZ%9vLf@v6rQ zk767mf#^(q2PY|Ww}-v;$}P4d?dMzPaf!8TgH5c~6Ftage4RaYx`y-U4xOXk%6 za}>*(zksCu7x*J;?T#`nR`D0(6&o1JDB4e>w%?|{eyvfn^$lRO`L8!> zzESfQMf4ku5Upv;qAh>3QOkzLE86O}8nvo7SkWfG)2PWhwkg{7cN?``%MwK!{z0RL z8?iM}B!AQh$p*|x6vrPo!m+mLh&KIVqo$3FMYPeE8a0x;t!TZkpk6qgJo}r-)3p9G zf`6Zo>>6h?(_q5VHAg%P?Z2Yh58zV?kx3esZb3n`lxAjF`i7cPq9`a8md4*sZ;Xf5 zM)3wkk{av^mU?IG=QrWtldEA4EyE~~$8@@Tg(I&ibiT?1isrrtPyv75r`Aqm> z8r(6kOx6ioU;D3NO6l6WB3-MB*GCb?!6in|RWr}w4N_Pvffhup^anYF0;kZ8YHrUU zihJIx_S=}M2Ca5Y^WFH`CA2(v{0Y?UQ}7kF_aEAJ6Je9%FB1OgxO5N~CQt5>_u%MSxx zE?Q5W6R!t@5ieA+b;;#L0MmAW+V51cL2BU&cKiGAF$cO}I4>``rX5t-rA3yD}6ygT(5>Zs^tr-{^585^-g;<8i9W>URJ=xEO5K2Jiz0X&%@<9BS? za!fm`ZpPCR*;h(l1( zZa7`Ag#9G0xbM>{4nmhJ8a^1|qJfjab=Y(}yr$!!QfUR~P`6?Nk?YJ?aN(-KkxtwQ z?cmZ-V2r;Tt+6I(3EYLN^-RkaJ-ysb&d~F~hfNQaszQyf9mQ?v8h_&niUqp+P>BK& zM29Nod&=2;;pY5BEEniNKYZyvrnZl^Q}{*=ehnTR?rWt>2KzKz04#6N?$hvIX&$F7 zaOFWCjb*?oYq%hRNW+M6+khBxNFmpve(0kbRz%Y;lu^L)YLlK#@^lm|g}qdgD?eJ+INVzF>9TXp@473!~cfd3i=djj54|nJ2W#^V-P76IyC);#{9L`uHXK5P~*21p>y# zFI^lP1$1oU?8VFH#wH%o9!ATF$rE*vtgsJwG|6jcAufre?;@UcNXorCqu_b!l>X6eOGj;)$_|^Hb<&bbNGT zCW(Fl(?;Klj5d8?1RrfcXTUac8Dp3tDz&qdm!6m!d*s55c46}3xlsTg9tCwH4__P= zondHaFOH0j<3nvDdB6d`7%_NZWHG#ur+aYF)I|kDL@+0N0+bB1za5+vJ`>3(u$_LTCpVB z^TicOd2xlpK@mC^H;Z6V#dQxI7{(=hRL#<`7b(-}rViTtENjESlGMv(eI+qFd$La( zoa@tO2m7>>b9231Wu>l4|Des4>`1Mv@mh=kh=8zvP+>Z`W*x3VB;qC!x_!=d3I^E` zTnmEPD?a?@cA!wB{rm+?U$Di+N8aY=`?SU31iGRBdV9}wX!zID)5E`Xn;mZCph*tq zU=n&I?ZP6@+v!DoH;=px;&H6eqNEbvS|p|{_7Q{TaVbsZ^{0C^4JMAJd}YZYQz381 z;XCtOMzfg2+%FkGlIS1o;|64@HwnJ#z&uN?BE6xp`-@(j%v{@4Seq&Hi{4x z4T?Ucg;JK6v0MUnTRW}|;eXIqvhcg8JtUcuoaq(FFQjAx+{WS+5gZt6M1vj)12@YL_YXlTxXEXE7sW zwUt3e4lyjF3BgV%7h883FJJjS#@&S&9HiZa_jmqARIogfFUNd@!3n z2RhudymA$qJaKtEH9SANwIfaV&P6>bclpH%Xe{)H0JT}jWv;axY!Q_Q^pg61< zaJ0Z?_;bsBD}?c?=vAZP5m<+AAVY14 zzC+zr4QF$+SU}Lcuz)0}hPlQaMIcsm5Mbh=hwT-o0wzh5L_DE?A$~%B5`zvQpqY9S z+#!KQ*U_-rTpiz&GcnyJVN2bP1kd$8Eu-w?G<|Udua4#zRai+ZaZ+|VOn(M5K7(If z;uFD8PvMbkV#EnpCUpj_;?;W1i|dhxA$hQuLe6@A0ryLl3ki>xG>`NTe1#W`#T9-o z^cCOOsgR1`W0%YpNty>ffpx;$IH7rvmCUZd+Ih1>Dti@lALEroM*eF_KS$01fB1k| zL6d=N_DcszxHijT7$exaV8X&n!q$miWK(VwZ(+167&K-?&6pd`gnpv2j|nHYN93UK zfI|ZSL4v#{%Bo4j!{Fyipt1sSdYBj`cX-g~AFQ=gkzs+0W%=O(IT%4LUL1xUU{tUu z4n8rn3AqWq8;4jEEFV&F2si7vv_GuqG^>#}R*b|cF>$B{FMlvZVVyua0T!%xgHXO0 zGKmdaNgSuf=skf_{2iagr^OQEz1joB6`ozmle8EuTG%$htJ83E9a5P2_|%^;$5J-q zQoGeNH<~J&<*nc~w1ibxs)SA&9R*~BLYs5F<*EAkJqT)SS@glu4*rykP5-Hb*b zI(6f->54EA%)10nZ{}VWV4l02BiRHJ0OuvAuk^;!4@KYu#CcSF)z}(&v;z4ETd7u- zX@qVxaX?<9FcdD+6@#21j~K5^;SXXK*bUP=&O&o(lG3dQ_-<%F;l|z`U*5pClm%*8 zqz%F~?_(Z>kQ~0{8ix(l)f&QnUe#g|fDj@Qy-r^mo9KUd@^Q$EKQMz&mvS(Jdm(em zU+2O-jfy^OKExe|6}eO5p+@8=WHdm%qo4FCrwh;MV9(Fmo`m3fvG=x7(2StiTo?EGjgwPiDatK!;7^uMwrcj)iF1!b#)gG=XkP(pq$%MWU zEpTHfxS62Zz1lF`G;E+kF9+dVXwIzAQqm`Khjbb+vy2HNF#Xf%h7#P~N@0PVnOrb1 z3Oqtg>O&40EturtnX+dH^%6}q4gH4;iUKl4c_s(kk9E>}hAcf=alnE;%|<*>@kX+~ zLJ{(7G84qpHTRlm;#V1cR6d4)Gv1N>YVceGa$4vcYntf2YO- zP9_nbhJ%|3aE7dqjW{|^D$SFO!L${y?5becqv>fh#D^H}6E3^t4cbATC#XGbSb$Ds=@;2SnH@x_ zK)Z4kD1mx;g$X(bq;)vCRX??sII&+tR^sFyk)<9bsh|{C8KG1^8fcQj5fmXig@`@N z4wH+5d_~&W&r+3HxKcn}TN{IIppFesCXvSBVN>i)Zm1rXXj3XES-|iGUyN7oiV038 zyRPl5&|=-{k&6LrWlz3`Lo9L$b(8qA9bOg=Eiq%8MK$u-da!DWy*ip&*h$edAH+;# zCHm=scwk^^FhuA)5#6(q9h~o(eOe&6mwxV@yEgYUti1rm3g(iyW%*h{)y|`8uONr; zh*hZpVQ}i;MXXCv+&aG^MYEF9C_%IL*3iRalU^Au&tfP&^wSf1Fe1fe>&eW16=X)8 z#0AkqH2!~qI)ztNplmPdnLW?kewruYy;EKm+Dq)uX|p56g*j~k@zfQgl8nc~Vj}C{ zqbb|zWT(F!zH~L!lN!LOl6}#~6-ka}SNV%5sfSbNP`*1_UbU+@OOYB%jiUsgZIZWx zNbseU0dPkcnv-xdWdq(K4jBpzg;W`!onm7zt6fPw4cHC}Xh_IAQ_lm^DnlF+`2N(( z0OA-(@RqtP_(*~B5r7960QDy{wc5K|V0Ll0 zd%E2va8%^Tu`LnG#!{@4EZJ_|x^3FVNs}~b)28j$G;NZmcAKVuoW@O>CjH-jU)?s1 zTc=Ll-+ME&y9a<20ZH-qRi8fG%$u1vZ{ECl^XARl;q#rZeS2XO{}=BL>ca7s^0s%) zY4~2v_QF=&|7he?E2A0EtF6J0x1MX|;ts1`Zx!R-ieC;Jwq0==t$44(5rs}QZ`Q3) zRGm^Q?lxz_#uRd`SJ7~Dp%rfxe$8&h+bWjl*wt3tbIxD%s=ieTkPhu2G|PSk6}AF9 zK-&U8YLuygSt5dfQgs`&!qzlZXq1Tp(UTLsV!PNO`ovDLOZ3aqvhU3~^R2i~(Wmp< z&dP#hP@RDtMm#KgDYU)7@x7pR6+|c8Sqk%36#D4C74Ok31Mr$-FG`xTP|vFSRonMk zcHCXHy!ptQ$0)YgUMt>S^~+Y(HbFLZ5qHUa(6ZuPa=51L0i)xA@SZ4~8-92=ph4{L z!0o`4ePaf?>J8tu%b2as3k`oh0yjXR@G&2EPgueAs&7m-eEeEobI3Ty95<%?>w$$D z1FA9jU2is=c>y%rj~_jLa>!7~deRbcEbh8yH;5~(coRr&*j3vKFrv<*`J?$`_yuRI zW=ovJEK}a}8ANk)&Jtuf-kNU&B)V4KQ{vj_Zv2*x!>_tHHf#|e?O*mN5!7AqH1V1AQQwDt5U(ZyfO_39W=1ygiVy+ zi2JFm=2s#}2r-IayF>*k`>@UQ@gmVqFu*Li*uAn@#F ziWS5Rn|0>NesGRw8)99n5pywC%r$k4sW$DpQ}yT33+Es>GPj%gLz;y&&|K9*$2Y(P zaLsBoAw`mnTk1^+i(G{nHGzTds99beOFS^6MiuD0JtwR?<%K5kQl(92-Xl=~EqmKR zn2>o0CPVa=b^T%)l-qWmfDuwCMQexKF|7{6@_RSs{0!4zV0#oRKNxUYdQm z_ApE5o%uY2Izy|3*^Ikq>6byhq`QVo&_N4mKj8!+*iI!;{N9?iIv$5LCia!)m3;Vvzj|*baGP5{vP6U5|!#Bh3bOy%lfHqu;O^m{2gz zT$9QvSB79V!hohK3JCo??MuKFU`oPgBFHf>jZD_G6dO98Ifp^+%@(QE95BOrw91hx zev4~Cfn68o*`oZNeVo)bZ#YnTm=|n%m|*A&Rt`2;BZb_Acnd3V0;q{92o`>|LM^bM z?0Qt9A+gT9@JP4XYEIY)qG_ky}X27n!j0N&0LFS@UA#vZN z=;gJJl|@qRcqD|-EleV0!<6bBs*{ma?4aCmSO*a^BoO;PmUs#KH%x^$;sLI?*l#lyECxNG)h!jjrMRQXdNH2yn{1y z)-?z8A>YS&H^cOPh6EutzH-u07nCvqNPluUbWy`fiuW-LTVZHM^$L_P_8Q)d53MdB z_C;YCGa&m8ASERW5)7Rt@nAawmweBvq63;Nlc*9(DK~N zHB*quBaJ)0@rjQ^lX78DuR5Wc!({~{%fnQrpGLwmAYp_Z$oI?kyr~mNZ|vCCL+s&- zOd28t&kJm;Q5IyO$bLzJb&H~cWw8`hd+s+ig>J+LWw*4`N}Rw_qkRG-sUdM0Th)Nc zStbc2%L43^O!1?z9?^cTO}$Pw^`H(IKlU%|&DdQ)&s@s_-yF)C`b za)Ut0iBj;4Ov)}a^=0fzSQaNWL&-X=sD-FlEyqxT$w&Jes=h}N32~vVgEt1x*XDxw zvKj>|ewZ52is(W~{7I8m9`aC=RKdlvrz1lv{-6lMdT^>xC|A6^8&vG7bFGp0?6Ba~ zYXw#JkrM?wc(4F&4%cec;WANz$}PDa5wgrkcj^w5+8|xN<)b;6;vXG^lII7&g{?RjW%0v3h_X@gbG8OebLk!{a=Gvq+SBGS{#@ot@`V}&~2 z+=%Z;py|R3tp=d1m=EY;Ufs9-jT_xj&m_sL>V!T%tu2w4&p4kJocs^El z227d>76)6d5{&Z++po0V_DfJklO?;?fGlFONSV4M0o^6vUjPt7g;P+~Z#4YH0OFwf z%)q1(AAhhAaVBvZU=*v~ooRcr z*2!i(l!6-eGZ9&^!nW8B!Qcn7rSj+iD>4HOoth0@)Ed7l3tjXZ3uMWRcV?lMSH><$ z>gZsKRsl1E-15}>us0but3rk$+sT9Rd(t`$$Q;y24V!@`>B1w8UyweEOg(ujj68W* zajx66mX2ib*mb{=;0R1Bz}doIz{wk%97I{3k=U46j_)8q&04@tAVrZ*;BWjA3OZO(>UmIsb?>OE7+Oit&xW_P`ro1biPcs5#CueE@jBK?JT`Z zumfQ(h9PZlyUM0hF>xSfTMDQle!Iyqfr?)N?k}6tELpId5GMF{Xkm51GQe4asO6MD z2kNlcU~+Yl_N>x5y@$*~$DeS*b#cNsHz z>jTSya=k81jn(JF0WMoPy1fKlYkcpcg173gpDJOCjOIb5zW{p}2fGx;0zp1BfOD>O z&61YuTcSM#Z-rCO9Lq{pc>&&$KmgJE>P`57Jd;g3?-^)rIz7O0#SvL~6qe(ARHe5X zz+qKSZI85eu(n{Z*$o(k^3;tXVSQ9>VOQ(e#=HvbDr|62=qb2jU!y}VA5uGQZ#aSk z988N4MJCbVycp2a>wF0A8&WF_=lmZ(1L4ED;~E z0QeL{rVO*O-8f}~k>RA~wCM+<^%?9qU^@GW*pmYbOTpAZrPG7yehoW)7M$I>2FUtg z#~B$kaD8HBb#fk!s0`baS+n_wk7W$2E#6^Oaq2DE*x1P|IG}4uEQzxKcxuiP)N@)x z`@u2>+MT*J`XRB#c`U^M+g}&g$jWqTq}y5TFx^BO2_gK~fW|#(%lIDMFe%4EOiqTT zG^tu}X2Wff_QHzxRcictP5R5`J&(qfVt3sw!A}^IG%uSBP4I6 z93rXRL%NiFTXy3b9KzB|NRqHidzICYX*DGvI8a%9$a%0F8wq^QXNRd|5>~vNHCFLS z9oi9;$AYl{Ca#OWjn%NvS{0&!K6K|QI8y-RO$V{+D2Z_^jkm%BVNnx8T!fBiol%pW@ z;gypX3%STrjI=kKmWRk ziM%MdU);QD6Ou=953rSV6v2;+WA2wE{De5Et2`(k(#h{hp$``|@E&omPQG8j{Eq zlbYY3z@AG;_?Y+s-Q@WM`0=FT4<@NkB&DAepVCcU5HITFkLYfcvyl7@d`TpRnc&YP z{!>9R@hZGH^!u84oqm;vGju|@w$Xw>(mx6fOSfU7z%c2EOuQrFUWO`%XjFou1zK8H zR*NRN8xzNFt=s|pbP^Z3-r2d5XlN$I%@uyVM{Tx<}`+59fRDVYN zSzYJPjd$qek0vbp^Qmn81_4Fh{Y&ciC&f=O{HMiV z)^heU;%9a8uZX{@lG+2xOpOomm7u|QeFoG!pI4<{L)T9{UwrnWskn_+;^&fe@z;{n zUr$nhT6{@U_VeOz=;YrNe@iF-w)h2|{5wgjUtHF(zbpQpPX3bkWu5$$8#?v(lhm(@ zf1rUso&dj+Ao^+o`-e&DABkVnP5v=7ol=Fro{;fRlGHy<+WfO5_0N;k*Tiq=#{WY6 zrcQo6DgT$^w>0p#6O{i-{Ei0xYjKhn$?qmDelI2G--s{kdjB?Aj{hzh*S{D4K^Ol= z@t<_^Ka2mOlmAtGLnr@-?mTgt9Zy_fy3X(GI^uucaDQBG4u6oOzPYSn|0|L1|1G|y z1^Exf|Iy{&7JsCY@ebwB@o*x9Om){QxSqp_6;4R3^Ifu3izFX1#hJKgp6*+aR#Di`le7$I%<=a@h4KB}}| zH#tDLgF1JJa(C$5os@fz&b^m%cj?^SNQ`f08VpL_#mV>>>sC+Lle5;=3n3n;coTkvTmW{+yx0C*YR2 zEmC?i z1DK-pM>2qEq?7hD$jAZCCaKHW3Xc(3kM!|u=80^zCo?7SzVu|R{0Pc+h72_7?K$$FlQWV-D|@E=z*(9?7O|Du z;Gt-?Pa!PQ{qN9n?l)-z>*3P?$GLIIHSOc)$6;0HycklQf_=QrI6B+fmL^`iFLmAbmAkGSy=NP*-e#Usa@=)WtKP0WgA#xB}RP3PT{A#Lv zNW0yU&lm7NBD4#%@koWB?9`hDsZuCkJcPp=UrRARJ zfRcL%e@K1CeTZYBOGm4^Xf=b@%kcl36zEA2GbA3k%%jjnVN>A2K^f>m1*;3YYq;zs zxl~=8L>CBqm#pQ!wy{{HQ^*2HcTYjyqYEeO7F+{-SE2bm4Qo+uIkcg5mjy!RaWzZg zGs);@DFg{K*T)lLC?;asD6~J`#uuGo%biB$pk&3pnX6t>OC1j1yf}b6>tmiTm2k3O;wU z0%l$BOjOd=VOIQq959TXu8W)Jp4KI`!;#_vl~V6X*85G`O}Qz9cT7l zR*QMHPZ^R4Z5rUW(eb;(_I*cWcf)Nov{p^P*Eh<&REU3LqbT0>3Fw!D$0ofN^EYjM zx0pZB&V2O0n(=?MQLd$o|HF-~co_s1?3Ds{`GZm>l9iP_s9!x< zkF3A!@V z+ZS>=h}F~mh4mA0o74T~x-oy(ryE~|V7lR(2=VRFiV6RsJXlLj8T}B%|%e{r?)P0A9 ziY@}*%H=k{Q5i@L{tBbG>7mVi6-kE%>a+f-jWR!7 zKtHij6sZQ3%-_KLDGkV%EjQDEH=WZBE5KD3H^1`^%_5ECtyJD`ZIs>VjQu8~xVg&v zB+8Y$;hj=*JVT!bUea(X3I+Pj{S1TLE!&{g+!iw9Tiw_2Q*=Lzzsv3${P#`%`#JnY zVEX4d^9B6H)t~!C0&ay8dU#z*NL@g1aEeqze*i=qvPe=jB>B1lF-}*)*WL~lQ7MXd zVW%SRR`QY)U0zbW&iHn|1RUSOH-O`tqj)D(g+om40(4_2L3oKE1L{{>AMVfj0V;X3 zNcslXp}ZS0LS%Rn5#Ymn6-D8wDaLoK2wFm6!bSHOYP56EeGs6c`w;y_ zNEiivPBtd#y}OZ|BCmnc`YVR7DFI;3Y6f@ zDLz5zmUMX&_qVkoP_~7_savQUhj>q0BL$W15ZT2M<*oQG_SPFY$5Xc?ieOoN5M8hG zj;pQ(5qlT=f^WMy`ChZl?@gU_xxrgtInYK#+NcfLeeHg%TcA6%8lCnVWvRO-<0ZKT zw*9NIu~$N&;z{jmo75Iwis4C?x1(|tHtQ^QlLKlkhX}sW0h}bLZj548vy7zs!}ud>?Z`5&t>RpJ#nm#H#p#W@Y-4pl~t3i*F?RPzmfYq@-)xO z1b>N;Y`?^Ex5f<8cUL?u?kiOLC-EmUYia4LD@tjq($XJWQA!jgrBd|b4hlR?#ZPWA z9Tav7V7whc%Jbl zVPd*{e7Z(^dhr3#UN|n{tEv}Ur56$K1re5LK}2&n$>F^HOXx;7cYg($c!y-9kwvgg zK0(!Tf1RmXqt*5qezwb)Lad+awXG}nYjQk$rufwhV?`&3xT5VkWd)C4F#(?edx3>MIla5Xg#zwo65R^RfR_u*=?xa z6@#rvZ=9(YPLeRQ_&|FB^)BP|fO21ZSqd=P7$YRmn#~x zlw>SkL`}BQ3j2vZcPG}aQIrr1dydRYvQh@j4i9%2hi6oa!-L#S#kQ)z2a|(Jb)i9TB;e^LgP*U3jDhYURH9P=(V>g^Q6=8t zxH&)e8vq^XhtV>`)DDp`gXcH#YxCgnxCDI>J!Igu67}5Akb%IkDm)g5S{|aY6#RyP z7-u34BSyfm9C1n^*P?#t1wbcMtO^LS2#-|zXdI#m!F(Qlhz(zi5B`94+$%kOOAs2q zeCC1%2-%lW^}#8gzFCx7bVT4G@i+~^Gn40MA0L?>HO6L)sp-kb#?FnNGY*f;Ab)ts zczkU3;^gI712v{cinC7|ljn_*;*-YsSn=GDG5W+5z1n6>PJ@84OH&hLqkxVT&rV!E zH&(o0Jc^dZ$ysA!?9$jQ(9TX8v*=PG8yiK-^Twso>9ZGsd*soviLu!yhk)w**ldyT zo}ZjHMvSSE>DjTfmnTN1jj7AiQ`SQ zFU+AcU>muNF-#Mc#@WfKC#T0QT%0v7PEMQ~1@O^PP&e}E#Hj2HLpwV$GIj}XP>)=~ zV{_bU5_o8YR8i7sJbrPMfYj#*{y#fAHd&;>oSiJrP9r;nu};rYm)!XAv6;~!V`O@4 zhS+j`dh!wtlNgIOs01`s1qazc*hKVYwf$Hk{-d_2z zih+Y~y=hp1;njQR8jxa+K+mfPbVJpm6K{Zcf9#wH!Il;2*(84JM^_?YjMPgVL-a{0 zeOB1pi(&9b5;nxYX44`?WIC3Hd@wrbMWI2*hzgeH!s*=6A%p%O96Z%);J?F%5A!d5 zxr!ilS~j2!Vjb0uiwnFGW)|?sHtpyn$8C)kjMErb-ohK<4Xc@37^2)9z9k7CIyq<< zP#Uy|!bP7JpZbu6zmE+V%|aflscr*FZur;`Hz2*oY@)S`r*f+p6us%*%MBDlMGZ6v z)uO;<{>B#8xh)aXGrK5d;_Rh5u)D_nL|erVbL~GBkq7!}6}*A=M2<2t6_oPKLcx!EZ1LBtx878Zs8Epub(5bJUMy$21CaM1V02 zIvj7%C>Y1N>Gh*Y(E}@rn#WKyjA_GLlDT8X12mC!2NPPCNSnt|epFJ>BrE_78_i>k z9MzW4gkTS77+y);;=vxQJ20UhHksqF6zqn<)tNS`SBGVw z(4ajCVLUlW8KKM*ucx^_XPLp*k7fWoLOTWUcos7Nnl~V!AK8_XHNuybF%EQ~f1u#> z;icwvtc2p_OXH`Huau}6C_B6(zJOD?Q%OZAZ=B<=h-i<7WkgH*npKV1c%lp;W+Gun z@lztkf)UU&PSgOK@3Eod226_+sED@^`QKt7drT{dN}d~!%`qY)K~nLdA+V`y4jVS? zBd{4}*wRpwFkY9v8Z!P5?$wC!Sj{GC-d+zNt%Q4>Y*!QEOIbj=RZ`zPQ5X z0*wfUNyFeA^X#hD$27N!&5$IiSFcbL ze@spTt@29U$43yS&s){N9=ytPdp7NG|HK&H4!LNLh# z?FP2&JC$AC@Ji*%v8z|nyf%*{sK%avZ-YRr>>$C!Lr>dlP9;oIDv5YP|E2f|`B4ly zg@9)2i7Dj?lRFddOIOlrS6A^pKiekT)HW;IGE5KehnC5>O+}mx-P6Zb=u9QyMe&=OanHyIijV!=QmMUw@FPzHogXe~cx$-Ol{d?@DB1>5QKM+}# z*{Cd8k}mqU5L#LgESy*hSESt2x^jjLm$V{TMLA_7<|JE5h+)~_6*=1qPan@}9BLU+ zGbSTlk<7ts2dPQ%1P=pZ)y@SWhrRt<$1ZZOz+SK7J*_qnpHGqgujkz@$}Nua%|K5qHSA$jENX?3Ud}Fon_et&haL%^5JM-lVX6AR8gZvDp{mMv4u0#X~jxBE5^4;^5TdI6$MatTx}zleeTz%)+D@#OQcOha+%6gJccj1t7RtOu(8p_xn%Oe%MSX24?TZ?u_g+_dHC z16xckg{0&izhHX>>OBYxm{;_Zr*F0)U^E${>h;HT06@=~GiS~a=W|d~$E8L&A)Dw; zTy4o@N;Z4hOWHwB@3;eoc5D33@j*B^cV`XmJlkU}ek5f4L*r*qF5|w*c*2pU<#u=B z7I$d;0!r|bk%|{p;O`wb0mi>eqyxnjyfpqC;2j`U0iGZK(HC10d2J<-?#W^%O?~5- bMZtwDs9;#Oe*)mk?zX!F{1Zn}DgXZf>>7UM diff --git a/docs/source/_build/doctrees/grogu.doctree b/docs/source/_build/doctrees/grogu.doctree index a099c16517f1827390b27583eb61faee6e86121f..a4f2663ad1e4762769b723bda285acf29fabef8f 100644 GIT binary patch literal 33107 zcmd^IdypJQdDorpt<#-OvW1N#YwTcaCGJkLFfYMX@o6;e55+G${+ht@gX2}R;aIh%_lCgq zkKgEBXw*Ej8XRW?9tV=aagH3XS@Z+of7O}cXvN;e87()VXM-?p322aAUmL zpQ#f305T^yGRyW=M6!JtblO*gmp_NUYw>pke@8)s{cMK&0&i~51B20Yhj=Pb5KY-V zSDPAaN-;tXmM5o_qUSQkjYVt2tT$_x!AX;pb-t-+ytzLbZ%WZP0WR-|N1aC1+DM5X zm&~r*ipL>jzQg3A6pu7b$8W`#bEIEnWJqbnm+|_AY_&bM9JWki`yJqvK#9+imB208 zt1%~jL2zb4aK=hzpfhvK^u*@pf{4zER(dY>a7=b(f2@8W#p)}z!D=DLQocFnt$1IC+oCqhaDf45c(JL50(Tl~FcI>8Me1 z8iH+N!Ew+>Tk%|)wTsj|CjbJU`1W2Dbl_Kns^$O|n-@=!(#HC4tQp>F2D*Sox)y3> z1ctNHrd3L~GedXEfdi`%Mt;WdUCHPRZ7VeAjcDGmn!(cb*248Cjxu#fy{B1kGRh)VUzO=(iA1nIB43sNcG1?(1@u%-kf-c2BxS!CB%+IGZ~8}-_$0HzuZ zC(s2}@4UrK$`_&AzQCajklSa57`rDKTk<>^;Q~PZ!~U3GC=NWSf@nC%S=x)$aF`j%K>)>{z7b@DJa%xN&6I3=5- zOyfv^;z|L%Qx6o%V%3alVfnmSi>xdH@xCBxHa$Nq2X)gAUmtnoeQ%O@-JZrvIYC>9 zm%d?qWgUn(Wg@VgEGJxIEfP@shwV!jDDvP3vBu(g$GT8K1k{RWWwC5FDt2aawsKk~ z5h;YGzh;GS6U*FS$^bxE%y<-SoNCK{Djq`#76!2gGi0Pl-Sdn$QkUjk7yD7>kt2eE z^t9lVLveX4E&{GrR#T60*sXXvZ6g7RsdkvjcmnQ7jj0k~GF1+}3b8E=TJCQU3QxVQ zg~*C(xndmS)|e^#zpxK+?FW`VPO61zWp1o0UL#D$j!8%}8Bfx*todfsPK*1?ZUj3T z4n>S)WPL9*S-Eb)@NO6voX|E*@(2v3*bLM6qDIxIIzA6Gjc8Gd@1+<@$^QGbvY%I( z4^rzx_Gif`V9mWpbk?a)n(lW~JHbWX9Ueti(o6Ah^CIl+YSiFyb#A7sUfpyWE%#vz z6^}z2!i!C-<$i=NnK?^eH@pQN0&T|Q4Nn$lWc31tI~QA<@wB)pujkN0Ol^=fnr(Qp zPPU_w2wq@?A>1au^d4j8=Yah{LSZrA(0(MohQ!c3xuW$8MWiIWBhjOaYxs+UaSeeY zVteuqKPT; zBSg}W!;5pjlJt!S${Xf8eJ`PKJeel$z$Mp;swWf7I#cy;_P=+e3Qk*!D))Od40qlB zD*->7n*yd%N@(TrOU3EG`&J(uta4&}= zuOu8-Qtrh5VQ_&xzUstGEslH9wJ?r_9k0dW2p?kzyElPf$%J?DvRzQJ{3nIYe3Sbs$G8cZQT$!NdKwxDsnHRD87UalKT|KO?)Kx zt+-MSQVNy4kL3O(T5>nZ-9*v-Rs3WWgZt~K$$3WZ+fiqax$mGJ1lfHj{cK={*yY<= zZVMIhq}16*r7M&u6fDUvp$T`0+y_vU(#E;szAsf3Z|MCVpZP&}@DRZ_GS0)-@qEKH zLN8OXvN&~PWd+$@JKV#2|e-TkgjI+P$A-;^S0$Hw4>7ad!vc^Ar2v+VgR@ai_S{?Wm7Ot6^i z*~L#XFTJ1C26!uf>0qu$3B;B`ty8fA`OvJe&-9R16!kvS2a~?>PlQi7-)T0-7a0$= z)vJZ~=>+3W=D%d?d?`h(`zabG?@R4UO$$dwJs4`cmvvU!Zo_TK)F}$SqNPCVz$>E(S(WNw!b2))5DPt)m z2FCIO=+UFGTuQ7Q3?L0|L`{QKU1-(}iejMw$Trr^b;~$+>&?r?%Wk=oZjfq7EZ#Op zQL;Ayg1bAr+-7Fc;UW+hRGWu{=S|)EKwFn zRc1|(s)x98aSuikT^Ge{?)z~Oh+KkN!^xK}X zA;{d={2G<$QW~2BT?H{^cBFj3?3_g3u4X6Mk&%BHky6KU0Qqz+8$s(I*YXJBy+i|3 z#Sc=6sbW?ts)~h9E4I1fk@(t*SFcBS@VK`V32hNt9qMlN0kbyzEN-4QNOHk}z6)W=GM`@TgRZJ7Wz!!a} z`VB31S{3{^=-R6a{t|0j@S*C{#c}+8-#BK~PTJ`tlae(P9j9ktX5ul>($!2TRNBnM z@1x)L%mhK^X5#P&N^~jB#I#xlDZ3zr6?Wkd(73Bz$UeKf;v8NJ+w##5JP7PO?*0OOl5mf)7$ChRx!bavQ zJ6P8~lM=he#ikRq?qXhEbLq_s1GZH(beIS@J#Pd7{wif}b#>lYUn-rLH|YPP+qO90 z7Lw9g()1SRn@QhzVC>16o1Khl+}5zNS~LHeVAH>$QIxp~tKznp_npI^V zQw{M-+Rkx0wVh*Zlxwoi+}E0*VRY@)1pOE^VG*gPX~O~ED4Ru6)&ff^D%M_A9L;BL z2hA}q@g4Y?z5{gSr9%06#Zf-F9h7&a7AF^tERTnvV7U0Aq|jY1j_xZNx;<@^xmyoK zw@pSs$!wFEOw^F1DOICrbM;UQY?tB4=p$(uMw_-*71T_R{-K zS*Tk%Q76*POIF=!EG?Xb|76w8jinQdRBA-^r5jGp8_lKbpPPKjRd`p*UumUM?^S(1 z=^GE8aXI&DI?X7#EVGqr>3u6fw10<75$Qfb!{i;VYqp*whF$`JI>+I>QuYi)VzB?hpk|0)A^vvlY6|lG*0lgKK7WbbS zOH$q@CWf!GE)!#JNIEbvxKM-~(-Wfm$m}|F?QvpWpogNH7y%_aG0-TIpq*4!^xUL} zRA6E_0tF`qpUE?|*LhxWPm(c04O4WK^>paw<*wf6rLKpfn->8kGcUt8oRsuXQO(L5 z^>7Nz3I|_sRz{f(>NF=qO-TSs3jdZKf8A`nz4zJpZ9Np-YzQcs+1Lp?E2&{+MbATe zNCoDDBhb%$=-XxkhK|N;U;)Fopt2eZ1ANN}JbdUGgiZxJ&8m1`$9z*6IN%FU5jXVOTuIn>v=Y`O zYIvLo;8DZ3MgJb721z=i20loksNuhf=`3paA1bk^ft3nJ4GA)ZqlRa(d0)qOgW8Z` ze|yM4eY(U8cMT|B;2^gO7Oqkx;cL}k>Bz>#HyoKamS1!B%+ks2hX(t{bzY!8G}xQ; zjR(x5TwE}j@Ihu8EWe|am=YIZ)VZB(P~=w_JNKqYOpue?0gP7UV&a)le(J1Arz)x)VRu6A6dE5a?c}Bxko7H9gRar|fY^y6VE#^wzam);5LaWFnB9 zf-JftbXW49UnrE17f1Qo42_;#xFtOloeL+RWL&rrT6+ma>_Sm=Mh~li8^=-5yJph4 zZkwAit|2;2%ur3zfRZAhswY6_RllzHDe?7CbW%b}qMDZ9)Dxj&*!_H3GUxY=aeEfJ5chD5mfYI?;~;2~maJj0`N>a*O5$4uxRePH>ds8Hz1%f=b1z)#y_V zyt-XubvZ<)YNzCYsqv^YK1s&+$vz$<(~;;=-Jq|O8&;wR|Nk;Tt^03~=yW|xs!QDR z6bTuRTfQp#>%587))^_Dh+Fs|h2oYU04R@Jo~9CuTUe=Z+>#(uIBq%K#96~(nZzDZ zi#WbO8@61X4O=>N?;N#o<)ixu{`!u%(D(22$uj1F6^=*8qOrHVgLv1ZPK@+%#CXy- zUL8O@){(E=S65%eU)G6Z(pq-Dsn*al69W2w;;o2xucTq}LBs*gtZ|1T<(aGR;B!7X zJN#a}D9e8ze?p7BHjJ1@SKY!-hY=;#Hihh*)Z>XHGJ3L@cWH!JsKXsKBd3egbmMkt zlCfHr2wB?}!M(0Hns@h&WN2wf`({&KEbalEEDs4{JF#2uJ=_1J7>AH_fbSa&# z{hc&I%IA`+70-P~(5>rp-(F57kORo4Fxjc<v}xZDq2KnjFG1$ocXA3Px|G^? zcS^ugN>@q&lyw$gCk=j=9q0+5BDSTo&KBsbGv8wQ zI2!!)^!X|JHUh`mLwSZl5|LNgdqKpp#;rKE-{It_^U;2$X;j)|7J|}uxcXrd6&|im zi~hRBn+jK@kRe>04(k)wCmC09^^Jk>o4+ALI`43Tv>3miytn)h{Q?eh9`ROEVKWL_%l7PYJSjc~#-4s$RI%MZ2(i5l7rT7V3xO)$ zQEjy1q2nzizGEqEDP`0;8cd1H$gkmQH>EzJ+^=*OsEh(5xAd9H$gs!^lk> z?}NGU{AyW7Jr&C1x*U&T;jJQj7C7Qmr9dDd){Gg2URhPe!^qNww0lv#Sx!U}jAa${ z_R(*cH5zc1)i7bF*;iG~z>2C~8A%h-YLB(i^sp%{-IGa41?khAz@aQ~7^_@$0&}It zKA`4Uj?1NjWj76C@EHq|Nn+ih`EL5o8qLb94kZZ2YTr*K06}U-wl=K>jaJ2!??Mu1 zOclR8GZfF!dV}ON)*y&Ar|yKUvwM!jGfk^m#wxL}N^y8LfJ{@(iV56lX~b99Vb}~# z9Y2oh07)`K+iNU({u)-)<kV9HEJW+| zsEMRmlnR^Ww-#z(6Vl^VoW?3-KA#;uUJ*A z;ETY(w(yzDW_)=Ct=7xh@^Z<2 zl6;o$k`aEn{Sv!wyV#ia4)n{uy*_ylLiWPU+;n~RY#*#n74S9#%wC{882Egpr3iX)(zMgPB(;u!#cH9L|mO7a^ zy1=|feB#oR5N|K%V|_fiuf|xh(!4&b945 literal 23697 zcmds935+CHdEVKXp8MF@V|Ybmk3JyoyXd;j_V|9{v2Na53uePAE|kIjU2+wr!S%}T}dD^}QKqZPm0 zs99dre4@GVvF1mbMK%_g7ec=gl&vP)hZ<$usa68ZYd+j$(^Nfjs^PM1KN4Ey$nm|T zxm~b_?U9F@MSF}DA}6X^+9x)(8TgwGqi&YB%uTCF6;W^KK5Ls0$bNpbDKHwYn}8oR z?O8TfN)WaO*|5ZOjSW>D&uXqU)=f5UMp59bHzLqS^b9(ckZ>Omj*z@+JeT_Y6A7yV}vHhC0QZYTps;=DOU+}8FSqWDF zE-ID5)fHeEt_aLlIK9h32;{Flv%j%!PY5o@Y_bo;?`-_u#(K@KG^#+1GiEMTB$)wB zP6$la>}!c$`#KP5zX^SQJb$byQ6!ZpcJX8?PxOt{j>ug+c zA?_?oPEn)hu%O-AHf!~&WeCz@%G!TaG~U%4jdyA^j)KekSi$is*0!d5LM=OTi;Y0T z0*5O_kqy;NCup)G0_k@cDN>s3uqa;!ya~t)5Wmh@ z@?IJ7OKjNLJ#BLsl!Z-R5uGpBcTSVN#?bdR4SypASwJEy3xhHQe_3v^Cq;ajrZ1(? zfwPDjLCW8a5oOzp#?Asad{ldq*Cr24t|j)mc0oG7WrSa+x@LW`MJHjI$)m zk3lYHiV)@)e^#`ACy(|~qFoR^N*%Y-Y=6`KS2EMzRs`DLf_D_y9OT#_M{f^?hrMZ zclY6s?&D*h68nyr^w&W855W8%#YUP92%_8qddL<+0{E=(yDk?$cNU*JAvoU~g358W zql8^mYg8krM8RJv5dvNzvR#t!%dHY-hhz|pm~;8OV_l$_t;uHM$CBxlZ9mW?xtISM zL&O-F!KM`W_mi1!x7kO}oNERy8 zM0^JQ_Sanm-p%w-9u}~hY*M$EfW*f>Y$h9pXR2}~1SaF9$S)JyqOj@SM<~4bSQ8O8 zm5T9GAoI1rjlo~qN5vQcp3jqR@{ucyS=sj@eF_l=w=@SBBe84-X5H2}%HDFHZ5d&s zMiw`^;2U0}wr&NXvEc^>1=29a%drLQs;<(_bs~g^>Vl8C0kV6JY9F<~M57g>TGewjU&t6QL zG#RE>F=P}NCiK1`uoS=Q{qk!t#V?3@d;>G0iWq&$U~1Dr{#{c8^+`6(n#b|?##UHO ziL@d`#3y{Ykl$^+_wKPww(^cC> zaRbxFXu()oGVXLRKa@j2IhHE=3CwBLIPRURVkLu2ZOwGNCfi@%X|f3vSljigQ+D7W*a%*uox0U@??+d>y0qi&^Hf09 z!@5;=Hk@+vB23xNX=!1(xkR|Mb6S|zHL~J=F3?I?IuXA--BUP^;ZQEzs#$;zlJ}t5fp+>EAsU_x+33e!4p)ZCXIs&eA1hr4s*;Bio|HRIDqp`vnaFNzIqEqP>)Kt|#22k{&{% zynswsQunLCrUxax5`DX0N1Kc8H}K1q()|}a>FSs;C#8;A-hg>6>e!~v6&*yZDRE1y zF*d{30ZTOaFeSl6LwN@_cy9Y#Y!*`H?43Pv`!v;_cNK7Q>6o~6H7(7d+kMOGQ$MGjUFpKW4R>nroT9}M_iQi1#l?~bIiNn%q z+SDj$rT1f>UuD;F+Fd-Dc;0{j43#l?i z-ESp)=_X{BitrRQ>Ru5B9a_N0Y`%>Zh8Xx~h0*7`#^@64rvSRA0sT`hn|wzk9qb=> zQB_Ku*w>T9^KH&YF%GqOTnI)iTw#mGGOMAlTip|p3)gX&(n8PClHelZ}YkVhi-WAk$282=TBUShEA!Zd&!r z{UIq5ts#xE5gDIqd{Cd~*jU00&CCo9qIt*3B=!k-Ff4-0;rb3u%hPjmavl4rlX>?? zdniY;GT{CsgqhOo8T6o{RW?Xq+QkVM%)mojn-{3+QLauiA<6?-F{H8M<-WQH!qnsy z;;(A*eclAj9#=Y`@w=UenBPnGunxI(l+{?tl~#7ejE`A?`T*EAKH4)V2 z56qX+4EoP5c>3LRlNY;a@*UM=N;3m|Ry!Ot)&opnw4kf7=F{$|^a}D#@$tc)_;^3n zR!b+6KZmkQD90bPj54f>SQ;@xY@!-fERt*+9wB33`NT3FI{9$PIL#3&Lhf(vrj0gS ztqkH&#Dg2Y9}$H1;pY%Rj2SU-QtQjl^rajC?I@qb3T(#n$7_7_`4qOBvDCxYGx%_( zu~HL^A|%-;7TW#ld8^XotIjn4KCJ=I5|XcDnkFCCBq7!`YG!1Y)~kMbOV@x2)40)C zH&%`E_;^fys!V7Ddix~M4Rf}cGla5Xv6k0dV}+`@ZdLIx?D-e40@&P@tIsYYPKa?d zmqOTm0AJ#`mW`i{XAFt1>P1Rne}($%Xr{9RzrDg<#o09)tBcO8x7$%v4^OYGa2El& zE@OGGf(1W3hRcq>VtFg1%3-vE1xzq}*^X*e!B>6`LNG!bfY3BK_)Z$_8hi_4rm#DZ$x#CVi7fAzKhx*Jdx71w*)V|F?v+L0&4iwCP!qachD~yE%B#e#!ptwq#%rTO-4L7j+|cvR_HS&SZGb zEkIG*``ur$d0Be%H!ejWMvn^D2Xb-E-Tpm`rg>3zEFJ7myKB6zhUSE(H-d&V!UNED z!CX8kbiMx3bp;Ch*epL3KwJ2J9osTdYJ->!tR{e0HNYRu1-Z^8nsvUdV$84W{x}>#ykp;9pWSz&O03etwj_t@vfxvb+&#AB-JbwS>A5!% zDz@Y`RI&?y?lH6z8}sfOP;>9WPily|_o5`T7teQ^>|yr-YC(|Q2kGha^hB?ZHQo19 zs{)zIWam8@)=-V(s73Y;@x}*GGtH6o?o9C>&oVBUJG2RdeJt$8H0s8Cw)nR5y<5hn z<&i4W>0LQ&x^JNL>{*GAeYBn7%EsOzB6An0k=Sy7CZQ)kmJ^{uLcV=ecG4O!h~=0L z0fRdWcuctg zRWhT%^m1MN75d?hW!0#Ir+@$Jd?ro@xR=9!tEyb|%tm?9s9>*q_1N~#hmO&67_FYX zIK=B}LBhCFF+oTMXZjp7eTcs(oNwwM=M|2J z`w=KqY?3>L1-Iox*+myEr{`o`AfbFxhf?6fDJDt2DVpxdhxTL^+TPAzAd&CT$vjdF z*CL)_G3q;rhVDs-2o9pbiz36$CnN93VWoS~oX>}nCmIPQC7KS)JRi@8(uYU{@Loic zUpF0(wY4OdEgl;IcRI#n3YC_4?BBU2V^Ox)Z#+hjMLhNg^hB?5JT|E;wicERp@s}& z{sYjv49mW7ufuP|&>fRF2=j4v@D58JvBDWb?7&$D56$@zNxYD$8pkP9ep|~3=~4|H zS;v|AlZ&Se95K8JRdAClL_8N7IISqdJFBuhe%A(XZY0esK33eZDsIpM2F^zrmIHtR z!P;Kjq24|;ZW`izK5lCPPjnd2i2@t~OHL~O-zP-bmp>sWO7~!v4VO&8HTP`s&HIow(${_QhP4 z%wO61cG8wzX0ByITVW^@rV3}CE>wf=l?2OnbI1HeFGZpIGIf&;Z{NJPghUXc#g3zK zo3`NzeN4g={-P+E8O}2wxzX+fnqo5FP2_kqA4*<$BB7+hlR@ZjOb~xo^xT*asgDpv zAdnw+#8Kl2XGvNxVvuzZjF}vX9-`r>#f^956OtD<-paA)N-$^hq2vihLP-k-<|Gje zJuAVS%ZJp5U<3kv3I?}kxy@@Y7krD!7@~?A9Xt7Sv1zixVl)z+C^wi(u|@RuEDSPGSv8P+$@?$TDg39^`By$lJt z^cv5wrc>dSTKbNiI~EhZiRxX>v2M6c9ApfAWW2zAvNFGpCBu?Zwbg6r6zII?q${Ck2a1 zTXs1al$jBXB@Bs2n^vL*)VmVodY=m@!rWV^o9tZRP#Z~bWb{?;OIbL+KNtDj@c*aK zw9D|noDU_>e@ZAR|GCG)@dxvv_2EMW68U~Ho&wy>*`iGWM&t+*a?330)3;c=hWLw8 ztWV_Akr&%M#*yf12JqQ@D0#AxP*SqluV!ZSsHk}=A5I@q5#W1~N`BpRWYczNUbe{Q zHgKn7WTQ}NiEM5_yGx5~2(pN5{(+w8HI8g15(B6$9mjSGVa)>Sb{W=uO)ed8A@z_` zQE^Ofi^p4U%R1*u%!sdU`+>8GmFIQ~%m3fWRKxk7r0jb-JyE1a zB8>dUM0dJ6JrS!K!bs6cUt#2b0jOAZ{vJK?Fp@v@9Y(f6rf5ISMAoDwkrJmheN)?{ zO34vk@+8?bS(Gd(Cxg?>-8}yq-xa};&SamsrCbScU^Fit%gh5PADjD{6d6u~Qc5-h zDD|5vYf4;4m-(BJ_d~qHI3UYjIyBb6X02b*od*$YMmDBN5+`0IB$J;MUvilYij;?t zX^ocRlz-~DwK&n4*IB=lG=TOZJ50xCLLA1&E?B7!HsGFxNCQ}i z{K%{}&+KF4aUHsSv9CdTGQg(GjewF2mXef$lW~biEyzaYpLq_04wAdVhu_gvJ+i)@P0Gh*Sbv`rTiki%e_GF{b? zr?5szA47CFFG{Sc#Y;*TuoGB{QzZhPY#?(duCd9$;zt9k)_Gj+IJ2LP+a`_?%GUy$ zwYD8_F+wh^P?Z)c)_P+T;HK%08l{myzm{1;$r!{@wg_v&kKYU~;g+~h$tIAFrV(u{ z{n!%Ed~;3uW625vKPcgSR7w@T-wQC5br82? zFz_~z;0GLWDzQQ!XV0|Ri2RZ&V*8Ou2sx#2ce9kJBpAyJaLd*TJY4rgH#5A1D?3SF z71cwlQSnPi{DoQvkPyQUt&-l7Ye^YrL!Ho}qps+yg!35YdX)%~SdPR)1;cI<#8B2A zahecXBde7*97<=2jGJ@x1RzK+#>Xa7lcBc?L_UFhqY#yNJu|>&X}m#myiF)#)u|y> z$K3~xvFW;1S9iFG7-S@=V<-T38XD}%_|6N8L*#XfY~NeD!V#etQuIbJS_v z5-xO^lC?+6Hqzwa_D3Mb8@S9_h}v5X+`8gfTSW>6VFq@gJ$~8oHlPF2E*vC7fpnmm zYEF%GmlB-v!AMzFw~&s{Lbj#~M(_?Wuq~{@?6MXx_r8=-FGg=0$qSAL;XXWFG&Dm#YoScFXgfkR#Bb|yn#-7#Jc}k+`9=FH{Q4W$OrctI$B-X* z=Lh0lI%?(8Q7o6X1YFuOaA^wd(iEGkAP+ZrT$x8Nuh13AclAJW2?h`$NU6PSB(+IA z+*oh0DcwHP5OjYLiloNo6!JxhkP^d}578-d%5f`U2J{j#FuWAjy|(UG69>Qx4`&LK z+GBc;{9>oc#$wpK>A@5fL|H+E8^r8LwnzB|!7c)avfrSv0Fh<~JR+%C(X zAc`$ed}e?6X!9&=B$D>Qt-@x7u>S&Z5=W$pCz4`rASWEJTy0dyWo(kyM;ZWxIx;E{ IYaGk}3*hM&ZU6uP diff --git a/docs/source/_build/doctrees/index.doctree b/docs/source/_build/doctrees/index.doctree index e84e375d7cc29c2260c1d4befcbc4919f761efea..1baf1f4640c837561e1a94f7369bbb3822e2a7d4 100644 GIT binary patch delta 197 zcmX?Wut=Y^fpzM1@r|tOm{ctFLyJ?3iuF?x^DRS$qs$C{9f*O3p~&nv%gDJ0-(m@1^HFQiADM@`Q>>z`H3mT`sqda R>7|=v86(&?hjTJA0RWe-6%YUb delta 221 zcmdlfdP9u0fpw}c_eNH8Mnwz#jQreG{glMK%+wrxm;Ca)oczR;V*QfT;*$8}{FKzq z35*f!^*qUCsqsmP#i{Y7MLAP?j5A6~3W}}t^?|~AxtS$Jse1WE>H5j(nYu}tc_7ie z#N1RpAm1Luuyqb|)t%DA2Qs86H7&I$H7_}}c#2;SPjWf3vFdPR%aM&uNmMFH)K5uG Z%goCx$;{6y)(4pd)QN1?=ItCzOaQl$QLO*~ diff --git a/docs/source/_build/html/.buildinfo b/docs/source/_build/html/.buildinfo index ed21c1e..7580025 100644 --- a/docs/source/_build/html/.buildinfo +++ b/docs/source/_build/html/.buildinfo @@ -1,4 +1,4 @@ # Sphinx build info version 1 -# This file records the configuration used when building these files. When it is not found, a full rebuild will be done. -config: 0f8ff1598a5b65221fc2ef8f29fdbeac +# This file hashes the configuration used when building these files. When it is not found, a full rebuild will be done. +config: d8f86c3651bae6bbb2cfa348920c9196 tags: 645f666f9bcd5a90fca523b33c5a78b7 diff --git a/docs/source/_build/html/_modules/grogu/useful.html b/docs/source/_build/html/_modules/grogu/useful.html index 15af420..8407880 100644 --- a/docs/source/_build/html/_modules/grogu/useful.html +++ b/docs/source/_build/html/_modules/grogu/useful.html @@ -13,7 +13,7 @@ - + @@ -89,15 +89,22 @@ # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE # SOFTWARE. +from itertools import permutations, product + import numpy as np from scipy.special import roots_legendre -from itertools import permutations, product + +# Pauli matrices +tau_x = np.array([[0, 1], [1, 0]]) +tau_y = np.array([[0, -1j], [1j, 0]]) +tau_z = np.array([[1, 0], [0, -1]]) +tau_0 = np.array([[1, 0], [0, 1]]) # define some useful functions
[docs] -def hsk(dh, k=(0, 0, 0)): +def hsk(H, ss, sc_off, k=(0, 0, 0)): """ One way to speed up Hk and Sk generation """ @@ -105,13 +112,12 @@ k.shape = (-1,) # are from the sisl source # this generates the list of phases - phases = np.exp(-1j * np.dot(np.dot(np.dot(dh.rcell, k), dh.cell), dh.sc.sc_off.T)) + phases = np.exp(-1j * 2 * np.pi * k @ sc_off.T) - HKU = np.einsum("abc,c->ab", dh.hup, phases) - HKD = np.einsum("abc,c->ab", dh.hdo, phases) - SK = np.einsum("abc,c->ab", dh.sov, phases) + HK = np.einsum("abc,a->bc", H, phases) + SK = np.einsum("abc,a->bc", ss, phases) - return HKU, HKD, SK
+ return HK, SK @@ -175,80 +181,190 @@ -
-[docs] -def make_atran(nauc, dirs="xyz", dist=1): +
+[docs] +def commutator(a, b): + "Shorthand for commutator" + return a @ b - b @ a
+ + + +
+[docs] +def tau_u(u): """ - Simple pair generator. Depending on the value of the dirs - argument sampling in 1,2 or 3 dimensions is generated. - If dirs argument does not contain either of x,y or z - a single pair is returend. + Pauli matrix in direction u. """ - if not (sum([d in dirs for d in "xyz"])): - return (0, 0, [1, 0, 0]) + u = u / np.linalg.norm(u) # u is force to be of unit length + return u[0] * tau_x + u[1] * tau_y + u[2] * tau_z
- dran = len(dirs) * [np.arange(-dist, dist + 1)] - mg = np.meshgrid(*dran) - dirsdict = dict() - for d in enumerate(dirs): - dirsdict[d[1]] = mg[d[0]].flatten() - for d in "xyz": - if not (d in dirs): - dirsdict[d] = 0 * dirsdict[dirs[0]] - ucran = np.array([dirsdict[d] for d in "xyz"]).T - atran = [] - for i, j in list(product(range(nauc), repeat=2)): - for u in ucran: - if (abs(i - j) + sum(abs(u))) > 0: - atran.append((i, j, list(u))) +# +
+[docs] +def crossM(u): + """ + Definition for the cross-product matrix. + Acting as a cross product with vector u. + """ + return np.array([[0, -u[2], u[1]], [u[2], 0, -u[0]], [-u[1], u[0], 0]])
+ - return atran
+
+[docs] +def RotM(theta, u, eps=1e-10): + """ + Definition of rotation matrix with angle theta around direction u. + """ + u = u / np.linalg.norm(u) + + M = ( + np.cos(theta) * np.eye(3) + + np.sin(theta) * crossM(u) + + (1 - np.cos(theta)) * np.outer(u, u) + ) + M[abs(M) < eps] = 0.0 # kill off small numbers + return M
-
-[docs] -def add(x, y): - """The sum of two numbers for testing. - This function adds to numbers together. I only created this for testing documentation and examples. +
+[docs] +def RotMa2b(a, b, eps=1e-10): + """ + Definition of rotation matrix rotating unit vector a to unit vector b. + Function returns array R such that R@a = b holds. + """ + v = np.cross(a, b) + c = a @ b + M = np.eye(3) + crossM(v) + crossM(v) @ crossM(v) / (1 + c) + M[abs(M) < eps] = 0.0 # kill off small numbers + return M
- Parameters - ---------- - x : float - First number - y : float - Second number added to `x` - Returns - ------- - sum : int - The sum of the inputs - See Also - -------- - numpy.add : Adds more than two numbers. +
+[docs] +def spin_tracer(M): + """ + Spin tracer utility. + This akes an operator with the orbital-spin sequence: + orbital 1 up, + orbital 1 down, + orbital 2 up, + orbital 2 down, + that is in the SPIN-BOX representation, + and extracts orbital dependent Pauli traces. + """ - Notes - ----- - We can create some latex notes here [1]_ : + M11 = M[0::2, 0::2] + M12 = M[0::2, 1::2] + M21 = M[1::2, 0::2] + M22 = M[1::2, 1::2] - .. math:: a + b = c + M_o = dict() + M_o["x"] = M12 + M21 + M_o["y"] = 1j * (M12 - M21) + M_o["z"] = M11 - M22 + M_o["c"] = M11 + M22 - References - ---------- - .. [1] https://numpydoc.readthedocs.io/en/latest/format.html + return M_o
- Examples - -------- - >>> add(1, 2) - 3 +
+[docs] +def parse_magnetic_entity(dh, atom=None, l=None, **kwargs): + """ + Function to define orbital indeces of a given magnetic entity. + dh: a sisl Hamiltonian object + atom: an integer or list of integers, defining atom (or atoms) in the unicell forming the magnetic entity + l: integer, defining the angular momentum channel + """ + # case where we deal with more than one atom defining the magnetic entity + if type(atom) == list: + dat = [] + for a in atom: + a_orb_idx = dh.geometry.a2o(a, all=True) + if ( + type(l) == int + ): # if specified we restrict to given l angular momentum channel inside each atom + a_orb_idx = a_orb_idx[[o.l == l for o in dh.geometry.atoms[a].orbitals]] + dat.append(a_orb_idx) + orbital_indeces = np.hstack(dat) + # case where we deal with a singel atom magnetic entity + elif type(atom) == int: + orbital_indeces = dh.geometry.a2o(atom, all=True) + if ( + type(l) == int + ): # if specified we restrict to given l angular momentum channel + orbital_indeces = orbital_indeces[ + [o.l == l for o in dh.geometry.atoms[atom].orbitals] + ] + + return orbital_indeces # numpy array containing integers labeling orbitals associated to a magnetic entity.
+ + + +
+[docs] +def blow_up_orbindx(orb_indices): + """ + Function to blow up orbital indeces to make SPIN BOX indices. """ - return x + y
+ return np.array([[2 * o, 2 * o + 1] for o in orb_indices]).flatten()
+ + + +
+[docs] +def calculate_exchange_tensor(pair): + o1, o2, o3 = pair["energies"] # o1=x, o2=y, o3=z + # dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]), + # dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), + # dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), + + J_ii = np.array([o2[-1], o3[0], o1[0]]) # xx, yy, zz + J_S = -0.5 * np.array([o3[1] + o3[2], o2[1] + o2[1], o1[1] + o1[2]]) # yz, zx, xy + D = 0.5 * np.array([o1[1] - o1[2], o2[2] - o2[1], o3[1] - o3[2]]) # x, y, z + return J_ii.sum() / 3, D, np.concatenate([J_ii[:2] - J_ii.sum() / 3, J_S]).flatten()
+ + + + diff --git a/docs/source/_build/html/_modules/index.html b/docs/source/_build/html/_modules/index.html index 85e47e5..a7e8580 100644 --- a/docs/source/_build/html/_modules/index.html +++ b/docs/source/_build/html/_modules/index.html @@ -13,7 +13,7 @@ - + diff --git a/docs/source/_build/html/_static/basic.css b/docs/source/_build/html/_static/basic.css index 93a4776..91cf78f 100644 --- a/docs/source/_build/html/_static/basic.css +++ b/docs/source/_build/html/_static/basic.css @@ -1,5 +1,12 @@ /* + * basic.css + * ~~~~~~~~~ + * * Sphinx stylesheet -- basic theme. + * + * :copyright: Copyright 2007-2024 by the Sphinx team, see AUTHORS. + * :license: BSD, see LICENSE for details. + * */ /* -- main layout ----------------------------------------------------------- */ @@ -108,11 +115,15 @@ img { /* -- search page ----------------------------------------------------------- */ ul.search { - margin-top: 10px; + margin: 10px 0 0 20px; + padding: 0; } ul.search li { - padding: 5px 0; + padding: 5px 0 5px 20px; + background-image: url(file.png); + background-repeat: no-repeat; + background-position: 0 7px; } ul.search li a { diff --git a/docs/source/_build/html/_static/doctools.js b/docs/source/_build/html/_static/doctools.js index 0398ebb..4d67807 100644 --- a/docs/source/_build/html/_static/doctools.js +++ b/docs/source/_build/html/_static/doctools.js @@ -1,5 +1,12 @@ /* + * doctools.js + * ~~~~~~~~~~~ + * * Base JavaScript utilities for all Sphinx HTML documentation. + * + * :copyright: Copyright 2007-2024 by the Sphinx team, see AUTHORS. + * :license: BSD, see LICENSE for details. + * */ "use strict"; diff --git a/docs/source/_build/html/_static/language_data.js b/docs/source/_build/html/_static/language_data.js index a5ea78e..434cd3d 100644 --- a/docs/source/_build/html/_static/language_data.js +++ b/docs/source/_build/html/_static/language_data.js @@ -1,6 +1,13 @@ /* + * language_data.js + * ~~~~~~~~~~~~~~~~ + * * This script contains the language-specific data used by searchtools.js, * namely the list of stopwords, stemmer, scorer and splitter. + * + * :copyright: Copyright 2007-2024 by the Sphinx team, see AUTHORS. + * :license: BSD, see LICENSE for details. + * */ var stopwords = ["a", "and", "are", "as", "at", "be", "but", "by", "for", "if", "in", "into", "is", "it", "near", "no", "not", "of", "on", "or", "such", "that", "the", "their", "then", "there", "these", "they", "this", "to", "was", "will", "with"]; diff --git a/docs/source/_build/html/_static/searchtools.js b/docs/source/_build/html/_static/searchtools.js index 2c774d1..b08d58c 100644 --- a/docs/source/_build/html/_static/searchtools.js +++ b/docs/source/_build/html/_static/searchtools.js @@ -1,5 +1,12 @@ /* + * searchtools.js + * ~~~~~~~~~~~~~~~~ + * * Sphinx JavaScript utilities for the full-text search. + * + * :copyright: Copyright 2007-2024 by the Sphinx team, see AUTHORS. + * :license: BSD, see LICENSE for details. + * */ "use strict"; @@ -13,7 +20,7 @@ if (typeof Scorer === "undefined") { // and returns the new score. /* score: result => { - const [docname, title, anchor, descr, score, filename, kind] = result + const [docname, title, anchor, descr, score, filename] = result return score }, */ @@ -40,14 +47,6 @@ if (typeof Scorer === "undefined") { }; } -// Global search result kind enum, used by themes to style search results. -class SearchResultKind { - static get index() { return "index"; } - static get object() { return "object"; } - static get text() { return "text"; } - static get title() { return "title"; } -} - const _removeChildren = (element) => { while (element && element.lastChild) element.removeChild(element.lastChild); }; @@ -65,13 +64,9 @@ const _displayItem = (item, searchTerms, highlightTerms) => { const showSearchSummary = DOCUMENTATION_OPTIONS.SHOW_SEARCH_SUMMARY; const contentRoot = document.documentElement.dataset.content_root; - const [docName, title, anchor, descr, score, _filename, kind] = item; + const [docName, title, anchor, descr, score, _filename] = item; let listItem = document.createElement("li"); - // Add a class representing the item's type: - // can be used by a theme's CSS selector for styling - // See SearchResultKind for the class names. - listItem.classList.add(`kind-${kind}`); let requestUrl; let linkUrl; if (docBuilder === "dirhtml") { @@ -120,10 +115,8 @@ const _finishSearch = (resultCount) => { "Your search did not match any documents. Please make sure that all words are spelled correctly and that you've selected enough categories." ); else - Search.status.innerText = Documentation.ngettext( - "Search finished, found one page matching the search query.", - "Search finished, found ${resultCount} pages matching the search query.", - resultCount, + Search.status.innerText = _( + "Search finished, found ${resultCount} page(s) matching the search query." ).replace('${resultCount}', resultCount); }; const _displayNextItem = ( @@ -145,7 +138,7 @@ const _displayNextItem = ( else _finishSearch(resultCount); }; // Helper function used by query() to order search results. -// Each input is an array of [docname, title, anchor, descr, score, filename, kind]. +// Each input is an array of [docname, title, anchor, descr, score, filename]. // Order the results by score (in opposite order of appearance, since the // `_displayNextItem` function uses pop() to retrieve items) and then alphabetically. const _orderResultsByScoreThenName = (a, b) => { @@ -255,7 +248,6 @@ const Search = { searchSummary.classList.add("search-summary"); searchSummary.innerText = ""; const searchList = document.createElement("ul"); - searchList.setAttribute("role", "list"); searchList.classList.add("search"); const out = document.getElementById("search-results"); @@ -326,7 +318,7 @@ const Search = { const indexEntries = Search._index.indexentries; // Collect multiple result groups to be sorted separately and then ordered. - // Each is an array of [docname, title, anchor, descr, score, filename, kind]. + // Each is an array of [docname, title, anchor, descr, score, filename]. const normalResults = []; const nonMainIndexResults = []; @@ -345,7 +337,6 @@ const Search = { null, score + boost, filenames[file], - SearchResultKind.title, ]); } } @@ -363,7 +354,6 @@ const Search = { null, score, filenames[file], - SearchResultKind.index, ]; if (isMain) { normalResults.push(result); @@ -485,7 +475,6 @@ const Search = { descr, score, filenames[match[0]], - SearchResultKind.object, ]); }; Object.keys(objects).forEach((prefix) => @@ -596,7 +585,6 @@ const Search = { null, score, filenames[file], - SearchResultKind.text, ]); } return results; diff --git a/docs/source/_build/html/genindex.html b/docs/source/_build/html/genindex.html index 112bda8..585a91c 100644 --- a/docs/source/_build/html/genindex.html +++ b/docs/source/_build/html/genindex.html @@ -13,7 +13,7 @@ - + @@ -71,16 +71,35 @@

Index

- A + B + | C | G | H | M + | P + | R + | S + | T
-

A

+

B

+
+ +

C

+ + +
@@ -132,8 +151,6 @@

M

+

P

+ + + +
+ +

R

+ + + +
+ +

S

+ + +
+ +

T

+ + +
+ diff --git a/docs/source/_build/html/grogu.html b/docs/source/_build/html/grogu.html index 23c92ec..ab23e48 100644 --- a/docs/source/_build/html/grogu.html +++ b/docs/source/_build/html/grogu.html @@ -14,9 +14,8 @@ - + - @@ -50,11 +49,19 @@
  • grogu.example module
  • grogu.jij module
  • grogu.useful module
  • Module contents
    • @@ -103,62 +110,46 @@

    grogu.useful module

    -
    -grogu.useful.add(x, y)[source]
    -

    The sum of two numbers for testing.

    -

    This function adds to numbers together. I only created this for testing documentation and examples.

    -
    -
    Parameters:
    -
      -

    • x (float) – First number

      • -

    • y (float) – Second number added to x

      • -
    -
    -
    Returns:
    -

    sum – The sum of the inputs

    -
    -
    Return type:
    -

    int

    -
    -
    -
    -

    See also

    -
    -
    numpy.add

    Adds more than two numbers.

    -
    -
    -
    -

    Notes

    -

    We can create some latex notes here [1] :

    -
    -\[a + b = c\]
    -

    References

    - -

    Examples

    -
    >>> add(1, 2)
-3
-
    -
    +
    +grogu.useful.RotM(theta, u, eps=1e-10)[source]
    +

    Definition of rotation matrix with angle theta around direction u.

    -
    -grogu.useful.hsk(dh, k=(0, 0, 0))[source]
    -

    One way to speed up Hk and Sk generation

    +
    +grogu.useful.RotMa2b(a, b, eps=1e-10)[source]
    +

    Definition of rotation matrix rotating unit vector a to unit vector b. +Function returns array R such that R@a = b holds.

    -
    -grogu.useful.make_atran(nauc, dirs='xyz', dist=1)[source]
    -

    Simple pair generator. Depending on the value of the dirs -argument sampling in 1,2 or 3 dimensions is generated. -If dirs argument does not contain either of x,y or z -a single pair is returend.

    +
    +grogu.useful.blow_up_orbindx(orb_indices)[source]
    +

    Function to blow up orbital indeces to make SPIN BOX indices.

    +
    + +
    +
    +grogu.useful.calculate_exchange_tensor(pair)[source]
    +
    + +
    +
    +grogu.useful.commutator(a, b)[source]
    +

    Shorthand for commutator

    +
    + +
    +
    +grogu.useful.crossM(u)[source]
    +

    Definition for the cross-product matrix. +Acting as a cross product with vector u.

    +
    + +
    +
    +grogu.useful.hsk(H, ss, sc_off, k=(0, 0, 0))[source]
    +

    One way to speed up Hk and Sk generation

    @@ -176,6 +167,39 @@ If dirs argument does not contain either of x,y or z a kset of a single k-pont at the origin is returend.

    +
    +
    +grogu.useful.parse_magnetic_entity(dh, atom=None, l=None, **kwargs)[source]
    +

    Function to define orbital indeces of a given magnetic entity. +dh: a sisl Hamiltonian object +atom: an integer or list of integers, defining atom (or atoms) in the unicell forming the magnetic entity +l: integer, defining the angular momentum channel

    +
    + +
    +
    +grogu.useful.print_atomic_indices(pair, magnetic_entities, dh)[source]
    +
    + +
    +
    +grogu.useful.spin_tracer(M)[source]
    +

    Spin tracer utility. +This akes an operator with the orbital-spin sequence: +orbital 1 up, +orbital 1 down, +orbital 2 up, +orbital 2 down, +that is in the SPIN-BOX representation, +and extracts orbital dependent Pauli traces.

    +
    + +
    +
    +grogu.useful.tau_u(u)[source]
    +

    Pauli matrix in direction u.

    +
    +

    Module contents

    diff --git a/docs/source/_build/html/index.html b/docs/source/_build/html/index.html index 6708b9c..1c2feba 100644 --- a/docs/source/_build/html/index.html +++ b/docs/source/_build/html/index.html @@ -14,9 +14,8 @@ - + - diff --git a/docs/source/_build/html/modules.html b/docs/source/_build/html/modules.html index 68a6a30..d7f7f4a 100644 --- a/docs/source/_build/html/modules.html +++ b/docs/source/_build/html/modules.html @@ -14,9 +14,8 @@ - + - @@ -84,11 +83,19 @@
  • grogu.example module
  • grogu.jij module
  • grogu.useful module
  • Module contents
    • diff --git a/docs/source/_build/html/objects.inv b/docs/source/_build/html/objects.inv index 45e7134..9bfc70c 100644 --- a/docs/source/_build/html/objects.inv +++ b/docs/source/_build/html/objects.inv @@ -2,6 +2,6 @@ # Project: Grogu # Version: # The remainder of this file is compressed using zlib. -xڕM -@ =E@z$Ķv~3MՊv^TV~` UVVujZ`nD_ltj&B)5>ڷ3h-084sMš`kk=2tV Ձ].9}L:-VTGr^BR&0TC|٠G{!İw$>N6:QYcϑ; - \ No newline at end of file +xڕMn F>Ruvt)jxb'fPuH);0 b~H",{dPaMv#[TYMIq2_f)J|:kֺ[Watv\a*0*`8܀ё55j2TA[P*اXM\XRGXl +B +k.ʐdz,h3tM/x_뽺@RP Z@iyt?Eq7 \ No newline at end of file diff --git a/docs/source/_build/html/py-modindex.html b/docs/source/_build/html/py-modindex.html index afa884e..f2c3645 100644 --- a/docs/source/_build/html/py-modindex.html +++ b/docs/source/_build/html/py-modindex.html @@ -13,7 +13,7 @@ - + diff --git a/docs/source/_build/html/search.html b/docs/source/_build/html/search.html index 61399b5..e8b7a3c 100644 --- a/docs/source/_build/html/search.html +++ b/docs/source/_build/html/search.html @@ -14,7 +14,7 @@ - + diff --git a/docs/source/_build/html/searchindex.js b/docs/source/_build/html/searchindex.js index 3d96987..3ca980e 100644 --- a/docs/source/_build/html/searchindex.js +++ b/docs/source/_build/html/searchindex.js @@ -1 +1 @@ -Search.setIndex({"alltitles": {"Contents:": [[1, null]], "Indices and tables": [[1, "indices-and-tables"]], "Module contents": [[0, "module-grogu"]], "Submodules": [[0, "submodules"]], "asd documentation": [[1, null]], "grogu package": [[0, null]], "grogu.example module": [[0, "module-grogu.example"]], "grogu.jij module": [[0, "module-grogu.jij"]], "grogu.useful module": [[0, "module-grogu.useful"]], "src": [[2, null]]}, "docnames": ["grogu", "index", "modules"], "envversion": {"sphinx": 64, "sphinx.domains.c": 3, "sphinx.domains.changeset": 1, "sphinx.domains.citation": 1, "sphinx.domains.cpp": 9, "sphinx.domains.index": 1, "sphinx.domains.javascript": 3, "sphinx.domains.math": 2, "sphinx.domains.python": 4, "sphinx.domains.rst": 2, "sphinx.domains.std": 2, "sphinx.ext.viewcode": 1}, "filenames": ["grogu.rst", "index.rst", "modules.rst"], "indexentries": {"add() (in module grogu.useful)": [[0, "grogu.useful.add", false]], "grogu": [[0, "module-grogu", false]], "grogu.example": [[0, "module-grogu.example", false]], "grogu.jij": [[0, "module-grogu.jij", false]], "grogu.useful": [[0, "module-grogu.useful", false]], "hsk() (in module grogu.useful)": [[0, "grogu.useful.hsk", false]], "make_atran() (in module grogu.useful)": [[0, "grogu.useful.make_atran", false]], "make_contour() (in module grogu.useful)": [[0, "grogu.useful.make_contour", false]], "make_kset() (in module grogu.useful)": [[0, "grogu.useful.make_kset", false]], "module": [[0, "module-grogu", false], [0, "module-grogu.example", false], [0, "module-grogu.jij", false], [0, "module-grogu.useful", false]]}, "objects": {"": [[0, 0, 0, "-", "grogu"]], "grogu": [[0, 0, 0, "-", "example"], [0, 0, 0, "-", "jij"], [0, 0, 0, "-", "useful"]], "grogu.useful": [[0, 1, 1, "", "add"], [0, 1, 1, "", "hsk"], [0, 1, 1, "", "make_atran"], [0, 1, 1, "", "make_contour"], [0, 1, 1, "", "make_kset"]]}, "objnames": {"0": ["py", "module", "Python module"], "1": ["py", "function", "Python function"]}, "objtypes": {"0": "py:module", "1": "py:function"}, "terms": {"0": 0, "1": 0, "150": 0, "2": 0, "20": 0, "3": 0, "42": 0, "A": 0, "If": 0, "One": 0, "The": 0, "ad": 0, "add": [0, 1, 2], "argument": 0, "b": 0, "c": 0, "can": 0, "contain": 0, "content": 2, "contour": 0, "creat": 0, "depend": 0, "detail": 1, "dh": 0, "dimens": 0, "dir": 0, "dist": 0, "document": 0, "doe": 0, "either": 0, "emax": 0, "emin": 0, "en": 0, "enum": 0, "exampl": 2, "first": 0, "float": 0, "format": 0, "function": 0, "gener": 0, "grid": 0, "grogu": [1, 2], "here": 0, "hk": 0, "hsk": [0, 2], "html": 0, "http": 0, "i": 0, "index": 1, "input": 0, "int": 0, "io": 0, "jij": 2, "k": 0, "kset": 0, "latest": 0, "latex": 0, "make_atran": [0, 2], "make_contour": [0, 2], "make_kset": [0, 2], "modul": [1, 2], "more": 0, "nauc": 0, "note": 0, "number": 0, "numk": 0, "numpi": 0, "numpydoc": 0, "onli": 0, "origin": 0, "p": 0, "packag": [1, 2], "page": 1, "pair": 0, "paramet": 0, "pont": 0, "readthedoc": 0, "refer": 0, "restructuredtext": 1, "returend": 0, "return": 0, "sampl": 0, "search": 1, "second": 0, "see": 1, "simpl": 0, "singl": 0, "sk": 0, "some": 0, "sophist": 0, "sourc": 0, "speed": 0, "src": 1, "submodul": 2, "sum": 0, "syntax": 1, "test": 0, "than": 0, "thi": 0, "togeth": 0, "two": 0, "type": 0, "up": 0, "us": [1, 2], "valu": 0, "wai": 0, "we": 0, "x": 0, "xyz": 0, "y": 0, "your": 1, "z": 0}, "titles": ["grogu package", "asd documentation", "src"], "titleterms": {"asd": 1, "content": [0, 1], "document": 1, "exampl": 0, "grogu": 0, "indic": 1, "jij": 0, "modul": 0, "packag": 0, "src": 2, "submodul": 0, "tabl": 1, "us": 0}}) +Search.setIndex({"alltitles": {"Contents:": [[1, null]], "Indices and tables": [[1, "indices-and-tables"]], "Module contents": [[0, "module-grogu"]], "Submodules": [[0, "submodules"]], "asd documentation": [[1, null]], "grogu package": [[0, null]], "grogu.example module": [[0, "module-grogu.example"]], "grogu.jij module": [[0, "module-grogu.jij"]], "grogu.useful module": [[0, "module-grogu.useful"]], "src": [[2, null]]}, "docnames": ["grogu", "index", "modules"], "envversion": {"sphinx": 62, "sphinx.domains.c": 3, "sphinx.domains.changeset": 1, "sphinx.domains.citation": 1, "sphinx.domains.cpp": 9, "sphinx.domains.index": 1, "sphinx.domains.javascript": 3, "sphinx.domains.math": 2, "sphinx.domains.python": 4, "sphinx.domains.rst": 2, "sphinx.domains.std": 2, "sphinx.ext.viewcode": 1}, "filenames": ["grogu.rst", "index.rst", "modules.rst"], "indexentries": {"blow_up_orbindx() (in module grogu.useful)": [[0, "grogu.useful.blow_up_orbindx", false]], "calculate_exchange_tensor() (in module grogu.useful)": [[0, "grogu.useful.calculate_exchange_tensor", false]], "commutator() (in module grogu.useful)": [[0, "grogu.useful.commutator", false]], "crossm() (in module grogu.useful)": [[0, "grogu.useful.crossM", false]], "grogu": [[0, "module-grogu", false]], "grogu.example": [[0, "module-grogu.example", false]], "grogu.jij": [[0, "module-grogu.jij", false]], "grogu.useful": [[0, "module-grogu.useful", false]], "hsk() (in module grogu.useful)": [[0, "grogu.useful.hsk", false]], "make_contour() (in module grogu.useful)": [[0, "grogu.useful.make_contour", false]], "make_kset() (in module grogu.useful)": [[0, "grogu.useful.make_kset", false]], "module": [[0, "module-grogu", false], [0, "module-grogu.example", false], [0, "module-grogu.jij", false], [0, "module-grogu.useful", false]], "parse_magnetic_entity() (in module grogu.useful)": [[0, "grogu.useful.parse_magnetic_entity", false]], "print_atomic_indices() (in module grogu.useful)": [[0, "grogu.useful.print_atomic_indices", false]], "rotm() (in module grogu.useful)": [[0, "grogu.useful.RotM", false]], "rotma2b() (in module grogu.useful)": [[0, "grogu.useful.RotMa2b", false]], "spin_tracer() (in module grogu.useful)": [[0, "grogu.useful.spin_tracer", false]], "tau_u() (in module grogu.useful)": [[0, "grogu.useful.tau_u", false]]}, "objects": {"": [[0, 0, 0, "-", "grogu"]], "grogu": [[0, 0, 0, "-", "example"], [0, 0, 0, "-", "jij"], [0, 0, 0, "-", "useful"]], "grogu.useful": [[0, 1, 1, "", "RotM"], [0, 1, 1, "", "RotMa2b"], [0, 1, 1, "", "blow_up_orbindx"], [0, 1, 1, "", "calculate_exchange_tensor"], [0, 1, 1, "", "commutator"], [0, 1, 1, "", "crossM"], [0, 1, 1, "", "hsk"], [0, 1, 1, "", "make_contour"], [0, 1, 1, "", "make_kset"], [0, 1, 1, "", "parse_magnetic_entity"], [0, 1, 1, "", "print_atomic_indices"], [0, 1, 1, "", "spin_tracer"], [0, 1, 1, "", "tau_u"]]}, "objnames": {"0": ["py", "module", "Python module"], "1": ["py", "function", "Python function"]}, "objtypes": {"0": "py:module", "1": "py:function"}, "terms": {"0": 0, "1": 0, "10": 0, "150": 0, "1e": 0, "2": 0, "20": 0, "3": 0, "42": 0, "A": 0, "If": 0, "One": 0, "act": 0, "add": 1, "ak": 0, "an": 0, "angl": 0, "angular": 0, "argument": 0, "around": 0, "arrai": 0, "atom": 0, "b": 0, "blow": 0, "blow_up_orbindx": [0, 2], "box": 0, "calculate_exchange_tensor": [0, 2], "channel": 0, "commut": [0, 2], "contain": 0, "content": 2, "contour": 0, "cross": 0, "crossm": [0, 2], "defin": 0, "definit": 0, "depend": 0, "detail": 1, "dh": 0, "dimens": 0, "dir": 0, "direct": 0, "doe": 0, "down": 0, "either": 0, "emax": 0, "emin": 0, "entiti": 0, "enum": 0, "ep": 0, "exampl": 2, "extract": 0, "form": 0, "function": 0, "gener": 0, "given": 0, "grid": 0, "grogu": [1, 2], "h": 0, "hamiltonian": 0, "hk": 0, "hold": 0, "hsk": [0, 2], "i": 0, "indec": 0, "index": 1, "indic": 0, "integ": 0, "jij": 2, "k": 0, "kset": 0, "kwarg": 0, "l": 0, "list": 0, "m": 0, "magnet": 0, "magnetic_ent": 0, "make": 0, "make_contour": [0, 2], "make_kset": [0, 2], "matrix": 0, "modul": [1, 2], "momentum": 0, "more": 0, "none": 0, "numk": 0, "object": 0, "oper": 0, "orb_indic": 0, "orbit": 0, "origin": 0, "p": 0, "packag": [1, 2], "page": 1, "pair": 0, "parse_magnetic_ent": [0, 2], "pauli": 0, "pont": 0, "print_atomic_indic": [0, 2], "product": 0, "r": 0, "represent": 0, "restructuredtext": 1, "returend": 0, "return": 0, "rotat": 0, "rotm": [0, 2], "rotma2b": [0, 2], "sampl": 0, "sc_off": 0, "search": 1, "see": 1, "sequenc": 0, "shorthand": 0, "simpl": 0, "singl": 0, "sisl": 0, "sk": 0, "sophist": 0, "sourc": 0, "speed": 0, "spin": 0, "spin_trac": [0, 2], "src": 1, "ss": 0, "submodul": 2, "syntax": 1, "tau_u": [0, 2], "theta": 0, "thi": 0, "trace": 0, "tracer": 0, "u": 0, "unicel": 0, "unit": 0, "up": 0, "us": [1, 2], "util": 0, "valu": 0, "vector": 0, "wai": 0, "x": 0, "xyz": 0, "y": 0, "your": 1, "z": 0}, "titles": ["grogu package", "asd documentation", "src"], "titleterms": {"asd": 1, "content": [0, 1], "document": 1, "exampl": 0, "grogu": 0, "indic": 1, "jij": 0, "modul": 0, "packag": 0, "src": 2, "submodul": 0, "tabl": 1, "us": 0}}) diff --git a/local_operator.ipynb b/local_operator.ipynb index 421af10..968674a 100644 --- a/local_operator.ipynb +++ b/local_operator.ipynb @@ -220,7 +220,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.7" + "version": "3.9.6" } }, "nbformat": 4, diff --git a/pyproject.toml b/pyproject.toml index 6be180a..187cb8e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -3,48 +3,52 @@ requires = ["hatchling"] build-backend = "hatchling.build" [project] -name = "grogu" -version = "0.0.1" +name = "grogu_magn" +version = "0.0.3" authors = [ - { name="Laszlo Oroszlany", email="laszlo.oroszlany@ttk.elte.hu" }, { name="Daniel Pozsar", email="danielpozsar@student.elte.hu" }, + { name="Laszlo Oroszlany", email="laszlo.oroszlany@ttk.elte.hu" }, ] description = "Relativistic magnetic interactions from non-orthogonal basis sets" readme = "README.md" kewords = [ "DFT", "physics", - "grogu" + "grogu", + "magnetic interactions", ] requires-python = ">=3.9" dependencies = [ - "numpy", + "numpy==1.24.4", "scipy", - "sisl", + "sisl==0.14.3", + "netcdf4==1.6.2", "openmpi", "mpi4py", - "argparse", ] + classifiers = [ "Development Status :: 3 - Alpha", "Programming Language :: Python", "Programming Language :: Python :: 3.9", - "Programming Language :: Python :: 3.12", "Programming Language :: Python :: 3 :: Only", "Topic :: Scientific/Engineering :: Physics", "License :: OSI Approved :: MIT License", "Operating System :: OS Independent", ] +[project.scripts] +grogu = "grogu:main" + [project.urls] -Homepage = "https://github.com/pypa/sampleproject" -Documentation = "https://readthedocs.org" -Repository = "https://github.com/me/spam.git" -Issues = "https://github.com/pypa/sampleproject/issues" +Homepage = "https://gitea.vo.elte.hu/et209d/grogu" +Documentation = "https://gitea.vo.elte.hu/et209d/grogu" +Repository = "https://gitea.vo.elte.hu/et209d/grogu" +Issues = "https://gitea.vo.elte.hu/et209d/grogu" [tool.pytest.ini_options] pythonpath = [ - "src/grogu/" + "src/grogu_magn/" ] addopts = [ "--import-mode=importlib", diff --git a/requirements-dev.txt b/requirements-dev.txt index 58f3049..1d5eb40 100644 --- a/requirements-dev.txt +++ b/requirements-dev.txt @@ -14,6 +14,5 @@ pytest pytest-randomly pytest-cov isort -sphinx tqdm matplotlib diff --git a/requirements.txt b/requirements.txt index 4bd7910..8020276 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,5 +1,4 @@ # Package dependencies -argparse numpy==1.24.4 scipy sisl==0.14.3 diff --git a/src/grogu/__init__.py b/src/grogu_magn/__init__.py similarity index 100% rename from src/grogu/__init__.py rename to src/grogu_magn/__init__.py diff --git a/src/grogu/example.py b/src/grogu_magn/example.py similarity index 100% rename from src/grogu/example.py rename to src/grogu_magn/example.py diff --git a/src/grogu_magn/grogu.py b/src/grogu_magn/grogu.py new file mode 100644 index 0000000..536734b --- /dev/null +++ b/src/grogu_magn/grogu.py @@ -0,0 +1,526 @@ +from useful import * + + +def main(): + import os + from sys import stdout + from timeit import default_timer as timer + + from tqdm import tqdm + + os.environ["OMP_NUM_THREADS"] = "1" # export OMP_NUM_THREADS=4 + os.environ["OPENBLAS_NUM_THREADS"] = "1" # export OPENBLAS_NUM_THREADS=4 + os.environ["MKL_NUM_THREADS"] = "1" # export MKL_NUM_THREADS=6 + os.environ["VECLIB_MAXIMUM_THREADS"] = "1" # export VECLIB_MAXIMUM_THREADS=4 + os.environ["NUMEXPR_NUM_THREADS"] = "1" # export NUMEXPR_NUM_THREADS=6 + + import warnings + + import numpy as np + import sisl + from mpi4py import MPI + from numpy.linalg import inv + + start_time = timer() + + # this cell mimicks an input file + fdf = sisl.get_sile("../../lat3_791/Fe3GeTe2.fdf") + # this information needs to be given at the input!! + scf_xcf_orientation = np.array([0, 0, 1]) # z + # list of reference directions for around which we calculate the derivatives + # o is the quantization axis, v and w are two axes perpendicular to it + # at this moment the user has to supply o,v,w on the input. + # we can have some default for this + ref_xcf_orientations = [ + dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]), + dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), + dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), + ] + + # human readable definition of magnetic entities + # magnetic_entities = [ + # dict(atom=0, ), + # dict(atom=1, ), + # dict(atom=2, ), + # dict(atom=3, l=2), + # dict(atom=4, l=2), + # dict(atom=5, l=2), + # ] + # pairs = [ + # dict(ai=3, aj=4, Ruc=np.array([0, 0, 0])), # isotropic should be -82 meV + # dict(ai=3, aj=5, Ruc=np.array([0, 0, 0])), # these should all be around -41.9 in the isotropic part + # dict(ai=4, aj=5, Ruc=np.array([0, 0, 0])), + # dict(ai=3, aj=0, Ruc=np.array([0, 0, 0])), + # dict(ai=3, aj=1, Ruc=np.array([0, 0, 0])), + # dict(ai=3, aj=2, Ruc=np.array([0, 0, 0])), + # ] + magnetic_entities = [ + dict(atom=3, l=2), + dict(atom=4, l=2), + dict(atom=5, l=2), + ] + + # pair information + pairs = [ + dict(ai=0, aj=1, Ruc=np.array([0, 0, 0])), # isotropic should be -82 meV + dict( + ai=0, aj=2, Ruc=np.array([0, 0, 0]) + ), # these should all be around -41.9 in the isotropic part + dict(ai=1, aj=2, Ruc=np.array([0, 0, 0])), + dict(ai=0, aj=1, Ruc=np.array([-1, 0, 0])), + dict(ai=0, aj=2, Ruc=np.array([-1, 0, 0])), + dict(ai=0, aj=1, Ruc=np.array([1, 0, 0])), + dict(ai=0, aj=2, Ruc=np.array([1, 0, 0])), + dict(ai=0, aj=1, Ruc=np.array([0, -1, 0])), + dict(ai=0, aj=2, Ruc=np.array([0, -1, 0])), + dict(ai=0, aj=1, Ruc=np.array([0, 1, 0])), + dict(ai=0, aj=2, Ruc=np.array([0, 1, 0])), + dict(ai=1, aj=2, Ruc=np.array([-1, 0, 0])), + ] + + # Brilloun zone sampling and Green function contour integral + kset = 20 + kdirs = "xy" + ebot = -30 + eset = 50 + esetp = 1000 + + # MPI parameters + comm = MPI.COMM_WORLD + size = comm.Get_size() + rank = comm.Get_rank() + root_node = 0 + if rank == root_node: + print("Number of nodes in the parallel cluster: ", size) + + simulation_parameters = dict( + path="Not yet specified.", + scf_xcf_orientation=scf_xcf_orientation, + ref_xcf_orientations=ref_xcf_orientations, + kset=kset, + kdirs=kdirs, + ebot=ebot, + eset=eset, + esetp=esetp, + parallel_size=size, + ) + + # digestion of the input + # read in hamiltonian + dh = fdf.read_hamiltonian() + try: + simulation_parameters["geom"] = fdf.read_geometry() + except: + print("Error reading geometry.") + + # unit cell index + uc_in_sc_idx = dh.lattice.sc_index([0, 0, 0]) + + setup_time = timer() + + NO = dh.no # shorthand for number of orbitals in the unit cell + + # preprocessing Hamiltonian and overlap matrix elements + h11 = dh.tocsr(dh.M11r) + h11 += dh.tocsr(dh.M11i) * 1.0j + h11 = h11.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") + + h22 = dh.tocsr(dh.M22r) + h22 += dh.tocsr(dh.M22i) * 1.0j + h22 = h22.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") + + h12 = dh.tocsr(dh.M12r) + h12 += dh.tocsr(dh.M12i) * 1.0j + h12 = h12.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") + + h21 = dh.tocsr(dh.M21r) + h21 += dh.tocsr(dh.M21i) * 1.0j + h21 = h21.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") + + sov = ( + dh.tocsr(dh.S_idx) + .toarray() + .reshape(NO, dh.n_s, NO) + .transpose(0, 2, 1) + .astype("complex128") + ) + + # Reorganization of Hamiltonian and overlap matrix elements to SPIN BOX representation + U = np.vstack( + [np.kron(np.eye(NO, dtype=int), [1, 0]), np.kron(np.eye(NO, dtype=int), [0, 1])] + ) + # This is the permutation that transforms ud1ud2 to u12d12 + # That is this transforms FROM SPIN BOX to ORBITAL BOX => U + # the inverse transformation is U.T u12d12 to ud1ud2 + # That is FROM ORBITAL BOX to SPIN BOX => U.T + + # From now on everything is in SPIN BOX!! + hh, ss = np.array( + [ + U.T + @ np.block([[h11[:, :, i], h12[:, :, i]], [h21[:, :, i], h22[:, :, i]]]) + @ U + for i in range(dh.lattice.nsc.prod()) + ] + ), np.array( + [ + U.T + @ np.block( + [[sov[:, :, i], sov[:, :, i] * 0], [sov[:, :, i] * 0, sov[:, :, i]]] + ) + @ U + for i in range(dh.lattice.nsc.prod()) + ] + ) + + # symmetrizing Hamiltonian and overlap matrix to make them hermitian + for i in range(dh.lattice.sc_off.shape[0]): + j = dh.lattice.sc_index(-dh.lattice.sc_off[i]) + h1, h1d = hh[i], hh[j] + hh[i], hh[j] = (h1 + h1d.T.conj()) / 2, (h1d + h1.T.conj()) / 2 + s1, s1d = ss[i], ss[j] + ss[i], ss[j] = (s1 + s1d.T.conj()) / 2, (s1d + s1.T.conj()) / 2 + + # identifying TRS and TRB parts of the Hamiltonian + TAUY = np.kron(np.eye(NO), tau_y) + hTR = np.array([TAUY @ hh[i].conj() @ TAUY for i in range(dh.lattice.nsc.prod())]) + hTRS = (hh + hTR) / 2 + hTRB = (hh - hTR) / 2 + + # extracting the exchange field + traced = [spin_tracer(hTRB[i]) for i in range(dh.lattice.nsc.prod())] # equation 77 + XCF = np.array( + [ + np.array([f["x"] for f in traced]), + np.array([f["y"] for f in traced]), + np.array([f["z"] for f in traced]), + ] + ) # equation 77 + + # Check if exchange field has scalar part + max_xcfs = abs(np.array(np.array([f["c"] for f in traced]))).max() + if max_xcfs > 1e-12: + warnings.warn( + f"Exchange field has non negligible scalar part. Largest value is {max_xcfs}" + ) + + H_and_XCF_time = timer() + + # for every site we have to store 3 Greens function (and the associated _tmp-s) in the 3 reference directions + for i, mag_ent in enumerate(magnetic_entities): + parsed = parse_magnetic_entity(dh, **mag_ent) # parse orbital indexes + magnetic_entities[i]["orbital_indeces"] = parsed + magnetic_entities[i]["spin_box_indeces"] = blow_up_orbindx( + parsed + ) # calculate spin box indexes + spin_box_shape = len( + mag_ent["spin_box_indeces"] + ) # calculate size for Greens function generation + + mag_ent["energies"] = ( + [] + ) # we will store the second order energy derivations here + + mag_ent["Gii"] = [] # Greens function + mag_ent["Gii_tmp"] = [] # Greens function for parallelization + mag_ent["Vu1"] = [ + list([]) for _ in range(len(ref_xcf_orientations)) + ] # These will be the perturbed potentials from eq. 100 + mag_ent["Vu2"] = [list([]) for _ in range(len(ref_xcf_orientations))] + for i in ref_xcf_orientations: + mag_ent["Gii"].append( + np.zeros((eset, spin_box_shape, spin_box_shape), dtype="complex128") + ) # Greens functions for every quantization axis + mag_ent["Gii_tmp"].append( + np.zeros((eset, spin_box_shape, spin_box_shape), dtype="complex128") + ) + + # for every site we have to store 2x3 Greens function (and the associated _tmp-s) + # in the 3 reference directions, because G_ij and G_ji are both needed + for pair in pairs: + spin_box_shape_i, spin_box_shape_j = len( + magnetic_entities[pair["ai"]]["spin_box_indeces"] + ), len( + magnetic_entities[pair["aj"]]["spin_box_indeces"] + ) # calculate size for Greens function generation + + pair["energies"] = [] # we will store the second order energy derivations here + + pair["Gij"] = [] # Greens function + pair["Gji"] = [] + pair["Gij_tmp"] = [] # Greens function for parallelization + pair["Gji_tmp"] = [] + + pair["Vij"] = [ + list([]) for _ in range(len(ref_xcf_orientations)) + ] # These will be the perturbed potentials from eq. 100 + pair["Vji"] = [list([]) for _ in range(len(ref_xcf_orientations))] + + for i in ref_xcf_orientations: + pair["Gij"].append( + np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype="complex128") + ) + pair["Gij_tmp"].append( + np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype="complex128") + ) # Greens functions for every quantization axis + pair["Gji"].append( + np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype="complex128") + ) + pair["Gji_tmp"].append( + np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype="complex128") + ) + + site_and_pair_dictionaries_time = timer() + + kset = make_kset(dirs=kdirs, NUMK=kset) # generate k space sampling + wkset = np.ones(len(kset)) / len(kset) # generate weights for k points + kpcs = np.array_split(kset, size) # split the k points based on MPI size + kpcs[root_node] = tqdm(kpcs[root_node], desc="k loop", file=stdout) + + k_set_time = timer() + + # this will contain all the data needed to calculate the energy variations upon rotation + hamiltonians = [] + + # iterate over the reference directions (quantization axes) + for i, orient in enumerate(ref_xcf_orientations): + # obtain rotated exchange field + R = RotMa2b(scf_xcf_orientation, orient["o"]) + rot_XCF = np.einsum("ij,jklm->iklm", R, XCF) + rot_H_XCF = sum( + [np.kron(rot_XCF[i], tau) for i, tau in enumerate([tau_x, tau_y, tau_z])] + ) + rot_H_XCF_uc = rot_H_XCF[uc_in_sc_idx] + + # obtain total Hamiltonian with the rotated exchange field + rot_H = hTRS + rot_H_XCF # equation 76 + + hamiltonians.append( + dict(orient=orient["o"], H=rot_H, rotations=[]) + ) # store orientation and rotated Hamiltonian + + for u in orient[ + "vw" + ]: # these are the infinitezimal rotations (for now) perpendicular to the quantization axis + Tu = np.kron(np.eye(NO, dtype=int), tau_u(u)) # section 2.H + + Vu1 = 1j / 2 * commutator(rot_H_XCF_uc, Tu) # equation 100 + Vu2 = 1 / 8 * commutator(commutator(Tu, rot_H_XCF_uc), Tu) # equation 100 + + for mag_ent in magnetic_entities: + mag_ent["Vu1"][i].append( + Vu1[:, mag_ent["spin_box_indeces"]][mag_ent["spin_box_indeces"], :] + ) # fill up the perturbed potentials (for now) based on the on-site projections + mag_ent["Vu2"][i].append( + Vu2[:, mag_ent["spin_box_indeces"]][mag_ent["spin_box_indeces"], :] + ) + + for pair in pairs: + ai = magnetic_entities[pair["ai"]][ + "spin_box_indeces" + ] # get the pair orbital sizes from the magnetic entities + aj = magnetic_entities[pair["aj"]]["spin_box_indeces"] + pair["Vij"][i].append( + Vu1[:, ai][aj, :] + ) # fill up the perturbed potentials (for now) based on the on-site projections + pair["Vji"][i].append(Vu1[:, aj][ai, :]) + + reference_rotations_time = timer() + + if rank == root_node: + print("Number of magnetic entities being calculated: ", len(magnetic_entities)) + print( + "We have to calculate the Greens function for three reference direction and we are going to calculate 15 energy integrals per site." + ) + print(f"The shape of the Hamiltonian and the Greens function is {NO}x{NO}.") + comm.Barrier() + # ---------------------------------------------------------------------- + + # make energy contour + # we are working in eV now ! + # and sisil shifts E_F to 0 ! + cont = make_contour(emin=ebot, enum=eset, p=esetp) + eran = cont.ze + + # ---------------------------------------------------------------------- + # sampling the integrand on the contour and the BZ + for k in kpcs[rank]: + wk = wkset[rank] # weight of k point in BZ integral + for i, hamiltonian_orientation in enumerate( + hamiltonians + ): # iterate over reference directions + # calculate Greens function + H = hamiltonian_orientation["H"] + HK, SK = hsk(H, ss, dh.sc_off, k) + Gk = inv(SK * eran.reshape(eset, 1, 1) - HK) + + # store the Greens function slice of the magnetic entities (for now) based on the on-site projections + for mag_ent in magnetic_entities: + mag_ent["Gii_tmp"][i] += ( + Gk[:, mag_ent["spin_box_indeces"]][..., mag_ent["spin_box_indeces"]] + * wk + ) + + for pair in pairs: + # add phase shift based on the cell difference + phase = np.exp(1j * 2 * np.pi * k @ pair["Ruc"].T) + + # get the pair orbital sizes from the magnetic entities + ai = magnetic_entities[pair["ai"]]["spin_box_indeces"] + aj = magnetic_entities[pair["aj"]]["spin_box_indeces"] + + # store the Greens function slice of the magnetic entities (for now) based on the on-site projections + pair["Gij_tmp"][i] += Gk[:, ai][..., aj] * phase * wk + pair["Gji_tmp"][i] += Gk[:, aj][..., ai] * phase * wk + + # summ reduce partial results of mpi nodes + for i in range(len(hamiltonians)): + for mag_ent in magnetic_entities: + comm.Reduce(mag_ent["Gii_tmp"][i], mag_ent["Gii"][i], root=root_node) + + for pair in pairs: + comm.Reduce(pair["Gij_tmp"][i], pair["Gij"][i], root=root_node) + comm.Reduce(pair["Gji_tmp"][i], pair["Gji"][i], root=root_node) + + green_function_inversion_time = timer() + + if rank == root_node: + # iterate over the magnetic entities + for tracker, mag_ent in enumerate(magnetic_entities): + # iterate over the quantization axes + for i, Gii in enumerate(mag_ent["Gii"]): + storage = [] + # iterate over the first and second order local perturbations + for Vu1, Vu2 in zip(mag_ent["Vu1"][i], mag_ent["Vu2"][i]): + # The Szunyogh-Lichtenstein formula + traced = np.trace( + (Vu2 @ Gii + 0.5 * Gii @ Vu1 @ Gii), axis1=1, axis2=2 + ) + # evaluation of the contour integral + storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we))) + + # fill up the magnetic entities dictionary with the energies + mag_ent["energies"].append(storage) + + # iterate over the pairs + for tracker, pair in enumerate(pairs): + # iterate over the quantization axes + for i, (Gij, Gji) in enumerate(zip(pair["Gij"], pair["Gji"])): + site_i = magnetic_entities[pair["ai"]] + site_j = magnetic_entities[pair["aj"]] + + storage = [] + # iterate over the first order local perturbations in all possible orientations for the two sites + for Vui in site_i["Vu1"][i]: + for Vuj in site_j["Vu1"][i]: + # The Szunyogh-Lichtenstein formula + traced = np.trace((Vui @ Gij @ Vuj @ Gji), axis1=1, axis2=2) + # evaluation of the contour integral + storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we))) + + # fill up the pairs dictionary with the energies + pairs[tracker]["energies"].append(storage) + + end_time = timer() + + print( + "############################### GROGU OUTPUT ###################################" + ) + print( + "================================================================================" + ) + print("Input file: ") + print(simulation_parameters["path"]) + print( + "Number of nodes in the parallel cluster: ", + simulation_parameters["parallel_size"], + ) + print( + "================================================================================" + ) + try: + print("Cell [Ang]: ") + print(simulation_parameters["geom"].cell) + except: + print("Geometry could not be read.") + print( + "================================================================================" + ) + print("DFT axis: ") + print(simulation_parameters["scf_xcf_orientation"]) + print("Quantization axis and perpendicular rotation directions:") + for ref in ref_xcf_orientations: + print(ref["o"], " --» ", ref["vw"]) + print( + "================================================================================" + ) + print("number of k points: ", simulation_parameters["kset"]) + print("k point directions: ", simulation_parameters["kdirs"]) + print( + "================================================================================" + ) + print("Parameters for the contour integral:") + print("Ebot: ", simulation_parameters["ebot"]) + print("Eset: ", simulation_parameters["eset"]) + print("Esetp: ", simulation_parameters["esetp"]) + print( + "================================================================================" + ) + print("Atomic informations: ") + print("") + print("") + print("Not yet specified.") + print("") + print("") + print( + "================================================================================" + ) + print("Exchange [meV]") + print( + "--------------------------------------------------------------------------------" + ) + print("Atom1 Atom2 [i j k] d [Ang]") + print( + "--------------------------------------------------------------------------------" + ) + for pair in pairs: + J_iso, J_S, D = calculate_exchange_tensor(pair) + J_iso = J_iso * sisl.unit_convert("eV", "meV") + J_S = J_S * sisl.unit_convert("eV", "meV") + D = D * sisl.unit_convert("eV", "meV") + + print(print_atomic_indices(pair, magnetic_entities, dh)) + print("Isotropic: ", J_iso) + print("DMI: ", D) + print("Symmetric-anisotropy: ", J_S) + print("") + + print( + "================================================================================" + ) + print("Runtime information: ") + print("Total runtime: ", end_time - start_time) + print( + "--------------------------------------------------------------------------------" + ) + print("Initial setup: ", setup_time - start_time) + print( + f"Hamiltonian conversion and XC field extraction: {H_and_XCF_time - setup_time:.3f} s" + ) + print( + f"Pair and site datastructure creatrions: {site_and_pair_dictionaries_time - H_and_XCF_time:.3f} s" + ) + print( + f"k set cration and distribution: {k_set_time - site_and_pair_dictionaries_time:.3f} s" + ) + print(f"Rotating XC potential: {reference_rotations_time - k_set_time:.3f} s") + print( + f"Greens function inversion: {green_function_inversion_time - reference_rotations_time:.3f} s" + ) + print( + f"Calculate energies and magnetic components: {end_time - green_function_inversion_time:.3f} s" + ) + + +if __name__ == "__main__": + main() diff --git a/src/grogu/jij.py b/src/grogu_magn/jij.py similarity index 99% rename from src/grogu/jij.py rename to src/grogu_magn/jij.py index 3de5d32..c91c289 100644 --- a/src/grogu/jij.py +++ b/src/grogu_magn/jij.py @@ -31,7 +31,7 @@ if __name__ == "__main__": from mpi4py import MPI from scipy.special import roots_legendre - from grogu.useful import hsk, make_atran, make_contour, make_kset + from grogu_magn.useful import hsk, make_atran, make_contour, make_kset start = timer() diff --git a/src/grogu/useful.py b/src/grogu_magn/useful.py similarity index 100% rename from src/grogu/useful.py rename to src/grogu_magn/useful.py diff --git a/test.ipynb b/test.ipynb index 7b1f469..07c01cd 100644 --- a/test.ipynb +++ b/test.ipynb @@ -9,7 +9,7 @@ "name": "stderr", "output_type": "stream", "text": [ - "[Daniels-Air:17787] shmem: mmap: an error occurred while determining whether or not /var/folders/yh/dx7xl94n3g52ts3td8qcxjcc0000gn/T//ompi.Daniels-Air.501/jf.0/2288320512/sm_segment.Daniels-Air.501.88650000.0 could be created.\n" + "[Daniels-Air:31350] shmem: mmap: an error occurred while determining whether or not /var/folders/yh/dx7xl94n3g52ts3td8qcxjcc0000gn/T//ompi.Daniels-Air.501/jf.0/1950875648/sm_segment.Daniels-Air.501.74480000.0 could be created.\n" ] } ], @@ -27,7 +27,7 @@ "\n", "import numpy as np\n", "import sisl\n", - "from grogu.useful import *\n", + "from grogu_magn.useful import *\n", "from mpi4py import MPI\n", "from numpy.linalg import inv\n", "import warnings\n", @@ -375,20 +375,7 @@ "Number of magnetic entities being calculated: 3\n", "We have to calculate the Greens function for three reference direction and we are going to calculate 15 energy integrals per site.\n", "The shape of the Hamiltonian and the Greens function is 84x84.\n", - "k loop: 27%|██▋ | 107/400 [00:29<01:20, 3.62it/s]\n" - ] - }, - { - "ename": "KeyboardInterrupt", - "evalue": "", - "output_type": "error", - "traceback": [ - "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", - "\u001b[0;31mKeyboardInterrupt\u001b[0m Traceback (most recent call last)", - "Cell \u001b[0;32mIn[7], line 26\u001b[0m\n\u001b[1;32m 24\u001b[0m H \u001b[38;5;241m=\u001b[39m hamiltonian_orientation[\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mH\u001b[39m\u001b[38;5;124m\"\u001b[39m]\n\u001b[1;32m 25\u001b[0m HK, SK \u001b[38;5;241m=\u001b[39m hsk(H, ss, dh\u001b[38;5;241m.\u001b[39msc_off, k)\n\u001b[0;32m---> 26\u001b[0m Gk \u001b[38;5;241m=\u001b[39m \u001b[43minv\u001b[49m\u001b[43m(\u001b[49m\u001b[43mSK\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43m \u001b[49m\u001b[43meran\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mreshape\u001b[49m\u001b[43m(\u001b[49m\u001b[43meset\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[43m \u001b[49m\u001b[43mHK\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 28\u001b[0m \u001b[38;5;66;03m# store the Greens function slice of the magnetic entities (for now) based on the on-site projections\u001b[39;00m\n\u001b[1;32m 29\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m mag_ent \u001b[38;5;129;01min\u001b[39;00m magnetic_entities:\n", - "File \u001b[0;32m<__array_function__ internals>:200\u001b[0m, in \u001b[0;36minv\u001b[0;34m(*args, **kwargs)\u001b[0m\n", - "File \u001b[0;32m~/Downloads/grogu/.venv/lib/python3.9/site-packages/numpy/linalg/linalg.py:538\u001b[0m, in \u001b[0;36minv\u001b[0;34m(a)\u001b[0m\n\u001b[1;32m 536\u001b[0m signature \u001b[38;5;241m=\u001b[39m \u001b[38;5;124m'\u001b[39m\u001b[38;5;124mD->D\u001b[39m\u001b[38;5;124m'\u001b[39m \u001b[38;5;28;01mif\u001b[39;00m isComplexType(t) \u001b[38;5;28;01melse\u001b[39;00m \u001b[38;5;124m'\u001b[39m\u001b[38;5;124md->d\u001b[39m\u001b[38;5;124m'\u001b[39m\n\u001b[1;32m 537\u001b[0m extobj \u001b[38;5;241m=\u001b[39m get_linalg_error_extobj(_raise_linalgerror_singular)\n\u001b[0;32m--> 538\u001b[0m ainv \u001b[38;5;241m=\u001b[39m \u001b[43m_umath_linalg\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43minv\u001b[49m\u001b[43m(\u001b[49m\u001b[43ma\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43msignature\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43msignature\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mextobj\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mextobj\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 539\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m wrap(ainv\u001b[38;5;241m.\u001b[39mastype(result_t, copy\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m))\n", - "\u001b[0;31mKeyboardInterrupt\u001b[0m: " + "k loop: 100%|██████████| 400/400 [01:46<00:00, 3.75it/s]\n" ] } ], @@ -453,9 +440,75 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 8, "metadata": {}, - "outputs": [], + "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", + "number of k points: 20\n", + "k point directions: xy\n", + "================================================================================\n", + "Parameters for the contour integral:\n", + "Ebot: -30\n", + "Eset: 50\n", + "Esetp: 10000\n", + "================================================================================\n", + "Atomic informations: \n", + "\n", + "\n", + "Not yet specified.\n", + "\n", + "\n", + "================================================================================\n", + "Exchange [meV]\n", + "--------------------------------------------------------------------------------\n", + "Atom1 Atom2 [i j k] d [Ang]\n", + "--------------------------------------------------------------------------------\n", + "[3]Fe(2) [4]Fe(2) [0 0 0] d [Ang] Not yet.\n", + "Isotropic: -60.89330922309355\n", + "DMI: [-5.95741551e+00 7.27737654e+00 6.90431275e-04 -8.11057566e-04\n", + " -5.49031203e-06]\n", + "Symmetric-anisotropy: [-9.32966923e-01 -8.92579299e-04 -2.04258659e-06]\n", + "\n", + "[3]Fe(2) [5]Fe(2) [0 0 0] d [Ang] Not yet.\n", + "Isotropic: -60.55651225519789\n", + "DMI: [-0.34565796 0.65919757 0.07106945 -6.23149871 -0.0424978 ]\n", + "Symmetric-anisotropy: [ 3.78506176e+00 -6.13838308e+00 3.59037036e-03]\n", + "\n", + "================================================================================\n", + "Runtime information: \n", + "Total runtime: 107.576339458\n", + "--------------------------------------------------------------------------------\n", + "Initial setup: 0.18185883300000016\n", + "Hamiltonian conversion and XC field extraction: 0.600 s\n", + "Pair and site datastructure creatrions: 0.010 s\n", + "k set cration and distribution: 0.016 s\n", + "Rotating XC potential: 0.230 s\n", + "Greens function inversion: 106.512 s\n", + "Calculate energies and magnetic components: 0.027 s\n" + ] + } + ], "source": [ "if rank == root_node:\n", " # iterate over the magnetic entities\n", @@ -565,7 +618,7 @@ }, { "cell_type": "code", - "execution_count": 25, + "execution_count": 9, "metadata": {}, "outputs": [], "source": [ @@ -575,16 +628,27 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "dh.geometry.plot(axes=\"xy\")" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 11, "metadata": {}, "outputs": [], "source": [ @@ -593,9 +657,27 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 12, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "xyz[-3:]: red, green, blue\n" + ] + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAABNoAAAHACAYAAAB0/gUQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABCsElEQVR4nO3df5hWdZ0//ucMyAwWM2rK8GMnf5b4E0iTMFu1Rln1YmU/5aLtAktqW9l+VLYfUiqWJqlpbEayaUbq9vG35ie9IGUjM1lNlNZMTRSElEH9qjOICjJzf//g47QToAOcmXtmeDyu61x5n/u87/v1PkznNfdzzn1ORalUKgUAAAAA2CqV5S4AAAAAAHoDQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEAB+pa7gO6otbU1zz//fAYMGJCKiopylwPQ45VKpaxatSpDhgxJZaW/8egzAMXSZzak1wAUq6O9RtC2Ec8//3zq6+vLXQZAr7N8+fL81V/9VbnLKDt9BqBz6DN/ptcAdI536zWCto0YMGBAkvU7r6ampszVAPR8zc3Nqa+vbzu+buv0GYBi6TMb0msAitXRXiNo24i3T62uqanRlAAK5Ksr6+kzAJ1Dn/kzvQagc7xbr3EBAwAAAAAogKANAAAAAAogaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgDQAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAAChA33IX0KusWJHMnp388Y/JgAHJCSckhx2WVFSUuzIAeoM330xuvDG59971veXww5NPfSqpri53ZQAA0C2tW5fccUcyd+76//7wh5N/+If1sU1nELQV5fLLkzPPTEqlpLLyz+uOOCK5/faktrac1QHQ0/32t8lxxyUvvpj0/X/t+6qrki99KbnrruRDHypvfQAA0M0sXpyMGZM888yff4X+8Y+TL385ueWW5Oiji39PXx0twq23Jv/7fyctLUlr6/qIdN269c/9+tfJ+PHlrQ+Anm3FiuSoo5KXX17/+H/2mZdeShoakhdeKF99AADQzbz+evLxjyfPPrv+8du/QpdK658bOzb5wx+Kf19B29YqlZJvfGPTXw9taVl/fuKiRV1aFgC9yKxZyWuvre8pf6mlJWlqWn92GwAAkCS54YZk+fKN/wrd2rp++bd/K/59BW1b67nnkv/+7/WB26b07bv+66MAsCVuumnjvyG8rbV1/TYAAECS9V8+rHyH1GvduvWXPy6aoG1rvf76u29TUZG88Ubn1wJA77R6dTHbAADANmL16vV/j34nb75Z/PsK2rZWfX3ynve88zZvvZXsv3/X1ANA7zNy5J+v3roxffuu3wYAAEiSDB/+zr9CV1Z2TlQjaNta/fsnJ5+c9Omz8ecrKpIddkhOOKFLywKgF/nCF/5884ONWbcu+fznu64eAADo5v75n9/5V+jW1uSLXyz+fQVtRfjmN5NhwzYM2/r0Wb/8x38k1dXlqQ2Anu+oo5LTTlv/3//z5jtvX3RiypTkiCO6vCwAAOiuhg1LLr54/X//z2u1VVSsX/7X/0r+8R+Lf19BWxFqa5Pf/CY566xkp53Wr6usXH+v2PvvT449trz1AdCzVVQkl1+e/PjHyb77/nn9/vsn116bfOc75asNAAC6qS9/OfnZz5KPfOTP63bfff3dRm+8cdNfTtwaFaXSO90uc9vU3Nyc2traNDU1paamZvMGt7YmTU3J9tsnVVWdUyBAD7NVx9VeaKv3R3Pz+vBtwIDiiwPogfSZDdknAO2tXr3+q6Q1Ne2/JNJRHT2uvsNl4dgilZXJjjuWuwoAejMfmAAAYLO8230si+KrowAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFCAsgZt9957b8aOHZshQ4akoqIit99+e7vnb7311hx99NF53/vel4qKiixatOhdX3P27NmpqKhot1RXV3fOBADo1vQZAACgK5U1aFu9enWGDx+emTNnbvL5ww47LBdddNFmvW5NTU1WrFjRtjz77LNFlAtAD6PPAAAAXalvOd/8mGOOyTHHHLPJ5ydMmJAkWbp06Wa9bkVFRQYNGrQ1pQHQC+gzAABAV+qV12h77bXXsuuuu6a+vj7HH398HnvssXfcfs2aNWlubm63AMCm6DMAAMDG9Lqgbe+9987VV1+dn/3sZ7nuuuvS2tqaQw89NH/60582OWb69Ompra1tW+rr67uwYgB6En0GAADYlF4XtI0ePToTJ07MiBEjcvjhh+fWW2/NLrvskn//93/f5JipU6emqampbVm+fHkXVgxAT6LPAAAAm1LWa7R1he222y4jR47M4sWLN7lNVVVVqqqqurAqAHoLfQYAAHhbrzuj7S+1tLTk0UcfzeDBg8tdCgC9kD4DAAC8raxntL322mvtzgBYsmRJFi1alJ122invf//78/LLL2fZsmV5/vnnkyRPPvlkkmTQoEFtd3ubOHFihg4dmunTpydJvvnNb+YjH/lI9tprr7z66qu55JJL8uyzz+aUU07p4tkBUG76DAAA0JXKGrQ99NBDOfLII9seT5kyJUkyadKkzJ49O3fccUcmT57c9vyJJ56YJJk2bVrOO++8JMmyZctSWfnnE/NeeeWVnHrqqWlsbMyOO+6Ygw46KPfff3/23XffLpgRAN2JPgMAAHSlilKpVCp3Ed1Nc3Nzamtr09TUlJqamnKXA9DjOa62Z38AFKu7H1fvvffeXHLJJVm4cGFWrFiR2267LePGjXvHMfPnz8+UKVPy2GOPpb6+PmeffXb+6Z/+qcPv2d33CUBP09Hjaq+/RhsAAEA5rV69OsOHD8/MmTM7tP2SJUty3HHH5cgjj8yiRYtyxhln5JRTTsncuXM7uVIAtlavv+soAABAOR1zzDE55phjOrz9rFmzsvvuu+fSSy9Nkuyzzz6577778t3vfjdjxozprDIBKIAz2gAAALqRBQsWpKGhod26MWPGZMGCBZscs2bNmjQ3N7dbAOh6gjYAAIBupLGxMXV1de3W1dXVpbm5OW+88cZGx0yfPj21tbVtS319fVeUCsBfELQBAAD0cFOnTk1TU1Pbsnz58nKXBLBNco02AACAbmTQoEFZuXJlu3UrV65MTU1N+vfvv9ExVVVVqaqq6oryAHgHzmgDAADoRkaPHp158+a1W3f33Xdn9OjRZaoIgI4StAEAAHSi1157LYsWLcqiRYuSJEuWLMmiRYuybNmyJOu/9jlx4sS27T/3uc/lmWeeyVe+8pU88cQT+cEPfpAbb7wxZ555ZjnKB2AzCNoAAAA60UMPPZSRI0dm5MiRSZIpU6Zk5MiROffcc5MkK1asaAvdkmT33XfPnXfembvvvjvDhw/PpZdemquuuipjxowpS/0AdJxrtAEAAHSiI444IqVSaZPPz549e6NjHnnkkU6sCoDO4Iw2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNgAAAAAogKANAAAAAAogaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgDQAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNgAAAAAogKANAAAAAAogaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgDQAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNgAAAAAogKANAAAAAAogaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgDQAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNgAAAAAogKANAAAAAAogaAMAAACAApQ1aLv33nszduzYDBkyJBUVFbn99tvbPX/rrbfm6KOPzvve975UVFRk0aJFHXrdm266KcOGDUt1dXUOOOCA3HXXXcUXD0C3p88AAABdqaxB2+rVqzN8+PDMnDlzk88fdthhueiiizr8mvfff39OOumknHzyyXnkkUcybty4jBs3Lr///e+LKhuAHkKfAQAAulJFqVQqlbuIJKmoqMhtt92WcePGbfDc0qVLs/vuu+eRRx7JiBEj3vF1xo8fn9WrV+fnP/9527qPfOQjGTFiRGbNmtWhWpqbm1NbW5umpqbU1NRszjQA2IjucFzVZwB6L8fVDdknAMXq6HG1112jbcGCBWloaGi3bsyYMVmwYMEmx6xZsybNzc3tFgDYGH0GAADYlF4XtDU2Nqaurq7durq6ujQ2Nm5yzPTp01NbW9u21NfXd3aZAPRQ+gwAALApvS5o2xJTp05NU1NT27J8+fJylwRAL6LPAADAtqFvuQso2qBBg7Jy5cp261auXJlBgwZtckxVVVWqqqo6uzQAegF9BgAA2JRed0bb6NGjM2/evHbr7r777owePbpMFQHQm+gzAADAppT1jLbXXnstixcvbnu8ZMmSLFq0KDvttFPe//735+WXX86yZcvy/PPPJ0mefPLJJOvPJnj7zIGJEydm6NChmT59epLk9NNPz+GHH55LL700xx13XK6//vo89NBD+eEPf9jFswOg3PQZAACgK5X1jLaHHnooI0eOzMiRI5MkU6ZMyciRI3PuuecmSe64446MHDkyxx13XJLkxBNPzMiRIzNr1qy211i2bFlWrFjR9vjQQw/NT3/60/zwhz/M8OHDc/PNN+f222/P/vvv34UzA6A70GcAAICuVFEqlUrlLqK7aW5uTm1tbZqamlJTU1PucgB6PMfV9uwPgGL1hOPqzJkzc8kll6SxsTHDhw/P5ZdfnkMOOWST28+YMSNXXHFFli1blp133jmf+tSnMn369FRXV3fo/XrCPgHoSTp6XO1112gDAADoTm644YZMmTIl06ZNy8MPP5zhw4dnzJgxeeGFFza6/U9/+tOcddZZmTZtWh5//PH86Ec/yg033JCvfe1rXVw5AJtL0AYAANCJLrvsspx66qmZPHly9t1338yaNSvbb799rr766o1uf//99+ejH/1oPv3pT2e33XbL0UcfnZNOOikPPvhgF1cOwOYStAEAAHSStWvXZuHChWloaGhbV1lZmYaGhixYsGCjYw499NAsXLiwLVh75plnctddd+XYY4/d5PusWbMmzc3N7RYAul5Z7zoKAADQm7300ktpaWlJXV1du/V1dXV54oknNjrm05/+dF566aUcdthhKZVKWbduXT73uc+941dHp0+fnm984xuF1g7A5nNGGwAAQDcyf/78XHjhhfnBD36Qhx9+OLfeemvuvPPOnH/++ZscM3Xq1DQ1NbUty5cv78KKAXibM9oAAAA6yc4775w+ffpk5cqV7davXLkygwYN2uiYc845JxMmTMgpp5ySJDnggAOyevXqfPazn83Xv/71VFZueL5EVVVVqqqqip8AAJvFGW0AAACdpF+/fjnooIMyb968tnWtra2ZN29eRo8evdExr7/++gZhWp8+fZIkpVKp84oFYKs5ow0AAKATTZkyJZMmTcrBBx+cQw45JDNmzMjq1aszefLkJMnEiRMzdOjQTJ8+PUkyduzYXHbZZRk5cmRGjRqVxYsX55xzzsnYsWPbAjcAuidBGwAAQCcaP358XnzxxZx77rlpbGzMiBEjMmfOnLYbJCxbtqzdGWxnn312KioqcvbZZ+e5557LLrvskrFjx+Zb3/pWuaYAQAdVlJx7vIHm5ubU1tamqakpNTU15S4HoMdzXG3P/gAoluPqhuwTgGJ19LjqGm0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFKGvQdu+992bs2LEZMmRIKioqcvvtt7d7vlQq5dxzz83gwYPTv3//NDQ05KmnnnrH1zzvvPNSUVHRbhk2bFgnzgKA7kyvAQAAukpZg7bVq1dn+PDhmTlz5kafv/jii/O9730vs2bNygMPPJD3vOc9GTNmTN588813fN399tsvK1asaFvuu+++zigfgB5ArwEAALpK33K++THHHJNjjjlmo8+VSqXMmDEjZ599do4//vgkyTXXXJO6urrcfvvtOfHEEzf5un379s2gQYM6pWYAeha9BgAA6Crd9hptS5YsSWNjYxoaGtrW1dbWZtSoUVmwYME7jn3qqacyZMiQ7LHHHvmHf/iHLFu27B23X7NmTZqbm9stAPR+XdVr9BkAANg2dNugrbGxMUlSV1fXbn1dXV3bcxszatSozJ49O3PmzMkVV1yRJUuW5GMf+1hWrVq1yTHTp09PbW1t21JfX1/MJADo1rqq1+gzAACwbei2QduWOuaYY3LCCSfkwAMPzJgxY3LXXXfl1VdfzY033rjJMVOnTk1TU1Pbsnz58i6sGICeZnN7jT4DAADbhm4btL193ZuVK1e2W79y5crNuibODjvskA9+8INZvHjxJrepqqpKTU1NuwWA3q+reo0+AwAA24ZuG7TtvvvuGTRoUObNm9e2rrm5OQ888EBGjx7d4dd57bXX8vTTT2fw4MGdUSYAPZheAwAAFKmsQdtrr72WRYsWZdGiRUnWX5R60aJFWbZsWSoqKnLGGWfkggsuyB133JFHH300EydOzJAhQzJu3Li21/jEJz6R73//+22Pv/SlL+VXv/pVli5dmvvvvz9/93d/lz59+uSkk07q4tkB0B3oNQAAQFfpW843f+ihh3LkkUe2PZ4yZUqSZNKkSZk9e3a+8pWvZPXq1fnsZz+bV199NYcddljmzJmT6urqtjFPP/10XnrppbbHf/rTn3LSSSfl//v//r/ssssuOeyww/Jf//Vf2WWXXbpuYgB0G3oNAADQVSpKpVKp3EV0N83NzamtrU1TU5Pr6AAUwHG1PfsDoFg94bg6c+bMXHLJJWlsbMzw4cNz+eWX55BDDtnk9q+++mq+/vWv59Zbb83LL7+cXXfdNTNmzMixxx7boffrCfsEoCfp6HG1rGe0AQAA9HY33HBDpkyZklmzZmXUqFGZMWNGxowZkyeffDIDBw7cYPu1a9fmqKOOysCBA3PzzTdn6NChefbZZ7PDDjt0ffEAbBZBGwAAQCe67LLLcuqpp2by5MlJklmzZuXOO+/M1VdfnbPOOmuD7a+++uq8/PLLuf/++7PddtslSXbbbbeuLBmALdRt7zoKAADQ061duzYLFy5MQ0ND27rKyso0NDRkwYIFGx1zxx13ZPTo0TnttNNSV1eX/fffPxdeeGFaWlo2+T5r1qxJc3NzuwWAridoAwAA6CQvvfRSWlpaUldX1259XV1dGhsbNzrmmWeeyc0335yWlpbcddddOeecc3LppZfmggsu2OT7TJ8+PbW1tW1LfX19ofMAoGMEbQAAAN1Ia2trBg4cmB/+8Ic56KCDMn78+Hz961/PrFmzNjlm6tSpaWpqaluWL1/ehRUD8DbXaAMAAOgkO++8c/r06ZOVK1e2W79y5coMGjRoo2MGDx6c7bbbLn369Glbt88++6SxsTFr165Nv379NhhTVVWVqqqqYosHYLM5ow0AAKCT9OvXLwcddFDmzZvXtq61tTXz5s3L6NGjNzrmox/9aBYvXpzW1ta2dX/84x8zePDgjYZsAHQfgjYAAIBONGXKlFx55ZX5yU9+kscffzyf//zns3r16ra7kE6cODFTp05t2/7zn/98Xn755Zx++un54x//mDvvvDMXXnhhTjvttHJNAYAO8tVRAACATjR+/Pi8+OKLOffcc9PY2JgRI0Zkzpw5bTdIWLZsWSor/3wORH19febOnZszzzwzBx54YIYOHZrTTz89X/3qV8s1BQA6qKJUKpXKXUR309zcnNra2jQ1NaWmpqbc5QD0eI6r7dkfAMVyXN2QfQJQrI4eV311FAAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACtC33AUA0L09+WTyk58kzz+fDBqUTJiQ7LdfuasCoLd45ZXkmmuS3/0uqapKxo5NxoxJ+vQpd2UAsPkEbQBsVGtrcvrpyfe/n/Ttm5RKSUVFctFFycknJ7NmrV8PAFvqlluSf/zHZM2aPwdrs2Yl+++fzJmTDB1a3voAYHP56igAG/Wtb60P2ZJk3bqkpWX9/ybJ1VcnX/96+WoDoOd74IFk/Pj1IVuptL7HvN1nnngiOfroPz8GgJ5ii4K2j3/84/nGN76xwfpXXnklH//4x7e6KADK6/XXk0su2fTzpVLyve8lTU2d8/76DEDv9+1vrz9TulTa8Ll165I//CG5666ur+ttkyZNyr333lu+AgDokbYoaJs/f36+//3vZ9y4cVm9enXb+rVr1+ZXv/pVYcUBUB733pusWvXO27z5ZnL33Z3z/voMQO/W0pL83//7zmes9e2b3HZb19X0l5qamtLQ0JAPfOADufDCC/Pcc8+VrxgAeowt/uroPffck8bGxnzkIx/J0qVLCywJgHJ7/fWObffGG51Xgz4D0Hu9fUmCd9La2rl95t3cfvvtee655/L5z38+N9xwQ3bbbbccc8wxufnmm/PWW2+VrzAAurUtDtoGDx6cX/3qVznggAPy4Q9/OPPnzy+wLADKaf/9i91uS+gzAL1XVVWy++7rvzr6Tjqzz3TELrvskilTpuR3v/tdHnjggey1116ZMGFChgwZkjPPPDNPPfVUeQsEoNvZoqCt4v91xKqqqvz0pz/N6aefnr/5m7/JD37wg0KLA6A8PvjB5PDD/3wHuL/Up09y0EHJyJGd8/76DEDv9y//8s7PV1Ymn/lM19TyblasWJG77747d999d/r06ZNjjz02jz76aPbdd99897vfLXd5AHQjfbdkUOkvrlh69tlnZ5999smkSZMKKQqA8rvqqmT06OTVV9tfQ6dv3+S9701+8pPOe299BqD3O+205Oc/T+bPX/810bf16bP+a6VXXJEMGVK28vLWW2/ljjvuyI9//OP84he/yIEHHpgzzjgjn/70p1NTU5Mkue222/KZz3wmZ555ZvkKBaBb2aKgbcmSJdlll13arfvkJz+ZYcOG5aGHHiqkMADKa6+9kocfTqZPT2bPXn+dnKqq5B//Mfna15I99ui899ZnAHq/fv3W31X03/4tufzy5E9/Wr/+r/86mTo1Oeqo8tY3ePDgtLa25qSTTsqDDz6YESNGbLDNkUcemR122KHLawOg+6oo/eVpA6S5uTm1tbVpampq+2sVwLZs3bqkuTkZMCDZbrvNH++42p79AdBeqZQ0Na0P37bffvPHd8Zx9dprr80JJ5yQ6urqQl6vq+k1AMXq6HF1i85oA2Db0rdvstNO5a4CgN6qoiLpbieGTZgwodwlANADbfFdRwEAAACAPxO0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABShr0Hbvvfdm7NixGTJkSCoqKnL77be3e75UKuXcc8/N4MGD079//zQ0NOSpp55619edOXNmdtttt1RXV2fUqFF58MEHO2kGAHR3eg0AANBVyhq0rV69OsOHD8/MmTM3+vzFF1+c733ve5k1a1YeeOCBvOc978mYMWPy5ptvbvI1b7jhhkyZMiXTpk3Lww8/nOHDh2fMmDF54YUXOmsaAHRjeg0AANBVKkqlUqncRSRJRUVFbrvttowbNy7J+jMMhgwZkn/913/Nl770pSRJU1NT6urqMnv27Jx44okbfZ1Ro0blwx/+cL7//e8nSVpbW1NfX59/+Zd/yVlnndWhWpqbm1NbW5umpqbU1NRs/eQAtnHd5bjaXXpNd9kfAL2F4+qG7BOAYnX0uNptr9G2ZMmSNDY2pqGhoW1dbW1tRo0alQULFmx0zNq1a7Nw4cJ2YyorK9PQ0LDJMQBsu/QaAACgSH3LXcCmNDY2Jknq6urara+rq2t77i+99NJLaWlp2eiYJ554YpPvtWbNmqxZs6btcXNz85aWDUAP0lW9Rp8BAIBtQ7c9o60rTZ8+PbW1tW1LfX19uUsCoBfRZwAAYNvQbYO2QYMGJUlWrlzZbv3KlSvbnvtLO++8c/r06bNZY5Jk6tSpaWpqaluWL1++ldUD0BN0Va/RZwAAYNvQbYO23XffPYMGDcq8efPa1jU3N+eBBx7I6NGjNzqmX79+Oeigg9qNaW1tzbx58zY5JkmqqqpSU1PTbgGg9+uqXqPPAADAtqGs12h77bXXsnjx4rbHS5YsyaJFi7LTTjvl/e9/f84444xccMEF+cAHPpDdd98955xzToYMGdJ2t7gk+cQnPpG/+7u/yxe/+MUkyZQpUzJp0qQcfPDBOeSQQzJjxoysXr06kydP7urpAdAN6DUAAEBXKWvQ9tBDD+XII49sezxlypQkyaRJkzJ79ux85StfyerVq/PZz342r776ag477LDMmTMn1dXVbWOefvrpvPTSS22Px48fnxdffDHnnntuGhsbM2LEiMyZM2eDi1YDsG3QawAAgK5SUSqVSuUuortpbm5ObW1tmpqafL0HoACOq+3ZHwDFclzdkH0CUKyOHle77TXaAAAAeouZM2dmt912S3V1dUaNGpUHH3ywQ+Ouv/76VFRUtLukAQDdl6ANAACgE91www2ZMmVKpk2blocffjjDhw/PmDFj8sILL7zjuKVLl+ZLX/pSPvaxj3VRpQBsLUEbAABAJ7rsssty6qmnZvLkydl3330za9asbL/99rn66qs3OaalpSX/8A//kG984xvZY489urBaALaGoA0AAKCTrF27NgsXLkxDQ0PbusrKyjQ0NGTBggWbHPfNb34zAwcOzMknn9yh91mzZk2am5vbLQB0PUEbAABAJ3nppZfS0tKywZ2p6+rq0tjYuNEx9913X370ox/lyiuv7PD7TJ8+PbW1tW1LfX39VtUNwJYRtAEAAHQTq1atyoQJE3LllVdm55137vC4qVOnpqmpqW1Zvnx5J1YJwKb0LXcBAAAAvdXOO++cPn36ZOXKle3Wr1y5MoMGDdpg+6effjpLly7N2LFj29a1trYmSfr27Zsnn3wye+655wbjqqqqUlVVVXD1AGwuZ7QBAAB0kn79+uWggw7KvHnz2ta1trZm3rx5GT169AbbDxs2LI8++mgWLVrUtvzt3/5tjjzyyCxatMhXQgG6OWe0AQAAdKIpU6Zk0qRJOfjgg3PIIYdkxowZWb16dSZPnpwkmThxYoYOHZrp06enuro6+++/f7vxO+ywQ5JssB6A7kfQBgAA0InGjx+fF198Meeee24aGxszYsSIzJkzp+0GCcuWLUtlpS8bAfQGFaVSqVTuIrqb5ubm1NbWpqmpKTU1NeUuB6DHc1xtz/4AKJbj6obsE4BidfS46s8mAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABSg2wdtq1atyhlnnJFdd901/fv3z6GHHprf/va3m9x+/vz5qaio2GBpbGzswqoB6En0GgAAoAh9y13AuznllFPy+9//Ptdee22GDBmS6667Lg0NDfnDH/6QoUOHbnLck08+mZqamrbHAwcO7IpyAeiB9BoAAKAI3fqMtjfeeCO33HJLLr744vz1X/919tprr5x33nnZa6+9csUVV7zj2IEDB2bQoEFtS2Vlt54qAGWi1wAAAEXp1p8I1q1bl5aWllRXV7db379//9x3333vOHbEiBEZPHhwjjrqqPzmN795x23XrFmT5ubmdgsA24au6DX6DAAAbBu6ddA2YMCAjB49Oueff36ef/75tLS05LrrrsuCBQuyYsWKjY4ZPHhwZs2alVtuuSW33HJL6uvrc8QRR+Thhx/e5PtMnz49tbW1bUt9fX1nTQmAbqYreo0+AwAA24aKUqlUKncR7+Tpp5/OZz7zmdx7773p06dPPvShD+WDH/xgFi5cmMcff7xDr3H44Yfn/e9/f6699tqNPr9mzZqsWbOm7XFzc3Pq6+vT1NTU7to7AGyZ5ubm1NbWdtvjamf3Gn0GoHN19z5TDvYJQLE6elzt1me0Jcmee+6ZX/3qV3nttdeyfPnyPPjgg3nrrbeyxx57dPg1DjnkkCxevHiTz1dVVaWmpqbdAsC2o7N7jT4DAADbhm4ftL3tPe95TwYPHpxXXnklc+fOzfHHH9/hsYsWLcrgwYM7sToAegO9BgAA2Bp9y13Au5k7d25KpVL23nvvLF68OF/+8pczbNiwTJ48OUkyderUPPfcc7nmmmuSJDNmzMjuu++e/fbbL2+++Wauuuqq/Od//md+8YtflHMaAHRjeg0AAFCEbh+0NTU1ZerUqfnTn/6UnXbaKZ/85CfzrW99K9ttt12SZMWKFVm2bFnb9mvXrs2//uu/5rnnnsv222+fAw88MPfcc0+OPPLIck0BgG5OrwEAAIrQ7W+GUA4uHApQLMfV9uwPgGI5rm7IPgEoVq+5GQIAAEBPN3PmzOy2226prq7OqFGj8uCDD25y2yuvvDIf+9jHsuOOO2bHHXdMQ0PDO24PQPchaAMAAOhEN9xwQ6ZMmZJp06bl4YcfzvDhwzNmzJi88MILG91+/vz5Oemkk/LLX/4yCxYsSH19fY4++ug899xzXVw5AJtL0AYAANCJLrvsspx66qmZPHly9t1338yaNSvbb799rr766o1u/x//8R/5whe+kBEjRmTYsGG56qqr0tramnnz5nVx5QBsLkEbAABAJ1m7dm0WLlyYhoaGtnWVlZVpaGjIggULOvQar7/+et56663stNNOm9xmzZo1aW5ubrcA0PUEbQAAAJ3kpZdeSktLS+rq6tqtr6urS2NjY4de46tf/WqGDBnSLqz7S9OnT09tbW3bUl9fv1V1A7BlBG0AAADd1Le//e1cf/31ue2221JdXb3J7aZOnZqmpqa2Zfny5V1YJQBv61vuAgAAAHqrnXfeOX369MnKlSvbrV+5cmUGDRr0jmO/853v5Nvf/nbuueeeHHjgge+4bVVVVaqqqra6XgC2jjPaAAAAOkm/fv1y0EEHtbuRwds3Nhg9evQmx1188cU5//zzM2fOnBx88MFdUSoABXBGGwAAQCeaMmVKJk2alIMPPjiHHHJIZsyYkdWrV2fy5MlJkokTJ2bo0KGZPn16kuSiiy7Kueeem5/+9KfZbbfd2q7l9t73vjfvfe97yzYPAN6doA0AAKATjR8/Pi+++GLOPffcNDY2ZsSIEZkzZ07bDRKWLVuWyso/f9noiiuuyNq1a/OpT32q3etMmzYt5513XleWDsBmErQBAAB0si9+8Yv54he/uNHn5s+f3+7x0qVLO78gADqFa7QBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUABBGwAAAAAUQNAGAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEABBG0AAAAAUIBuH7StWrUqZ5xxRnbdddf0798/hx56aH7729++45j58+fnQx/6UKqqqrLXXntl9uzZXVMsAD2SXgMAABSh2wdtp5xySu6+++5ce+21efTRR3P00UenoaEhzz333Ea3X7JkSY477rgceeSRWbRoUc4444yccsopmTt3bhdXDkBPodcAAABFqCiVSqVyF7Epb7zxRgYMGJCf/exnOe6449rWH3TQQTnmmGNywQUXbDDmq1/9au688878/ve/b1t34okn5tVXX82cOXM69L7Nzc2pra1NU1NTampqtn4iANu47nxcLUev6c77A6AnclzdkH0CUKyOHle79Rlt69atS0tLS6qrq9ut79+/f+67776NjlmwYEEaGhrarRszZkwWLFiwyfdZs2ZNmpub2y0AbBu6otfoMwAAsG3o1kHbgAEDMnr06Jx//vl5/vnn09LSkuuuuy4LFizIihUrNjqmsbExdXV17dbV1dWlubk5b7zxxkbHTJ8+PbW1tW1LfX194XMBoHvqil6jzwAAwLahWwdtSXLttdemVCpl6NChqaqqyve+972cdNJJqawsrvSpU6emqampbVm+fHlhrw1A99fZvUafAQCAbUPfchfwbvbcc8/86le/yurVq9Pc3JzBgwdn/Pjx2WOPPTa6/aBBg7Jy5cp261auXJmampr0799/o2OqqqpSVVVVeO0A9Ayd3Wv0GQAA2DZ0+zPa3vae97wngwcPziuvvJK5c+fm+OOP3+h2o0ePzrx589qtu/vuuzN69OiuKBOAHkyvAQAAtka3D9rmzp2bOXPmZMmSJbn77rtz5JFHZtiwYZk8eXKS9V/HmThxYtv2n/vc5/LMM8/kK1/5Sp544on84Ac/yI033pgzzzyzXFMAoJvTawAAgCJ0+6Ctqakpp512WoYNG5aJEyfmsMMOy9y5c7PddtslSVasWJFly5a1bb/77rvnzjvvzN13353hw4fn0ksvzVVXXZUxY8aUawoAdHN6DQAAUISKUqlUKncR3U1zc3Nqa2vT1NSUmpqaDo9bsWpFZi+anT++/McM6DcgJ+x7Qg57/2GpqKjoxGoBur8tPa72Vlu6P95c92ZufOzG3PvsvalIRQ7f7fB8at9PpbpvdSdWC9D96TMb2pJ9Umptzb3/9/LcvOBHeW3dGxm24175p09fkrrd9+/kagG6v44eV7v9GW09xeUPXJ7679bn7F+enev++7pc8dAV+evZf52PX/PxNL3ZVO7yAOjhfvvcb/P+774/k26flJ/87ieZ/bvZmXDbhOw2Y7c8vOLhcpcHwLuYOXNmdtttt1RXV2fUqFF58MEH33H7m266KcOGDUt1dXUOOOCA3HXXXZ1a3ysrluSvp+yQIxadkVn9Hs112y/O19bOyV/9+IDMuuzTnfreAL2JoK0Atz5+a/73nP+dllJLWkutWde6Luta1yVJfv3srzP+5vFlrhCAnmzFqhU56tqj8vIbLydJuz7z0usvpeGahryw+oVylgjAO7jhhhsyZcqUTJs2LQ8//HCGDx+eMWPG5IUXNn7svv/++3PSSSfl5JNPziOPPJJx48Zl3Lhx+f3vf98p9ZVaW/PJiz6UBbWrkiTr+qxfWivX/+/nV/2f/Py6czvlvQF6G0HbViqVSvnGr76Rimz866EtpZbMfXpuFjUu6trCAOg1Zj00K6+tfS0tpZYNnmsptaRpTVOueviqMlQGQEdcdtllOfXUUzN58uTsu+++mTVrVrbffvtcffXVG93+3/7t3/I3f/M3+fKXv5x99tkn559/fj70oQ/l+9//fqfU99t7fpJf7vhqWjbx6bCyNblg4WWd8t4AvY2gbSs9t+q5/PfK/04pm77UXd+Kvrn9idu7rigAepWb/nDTRkO2t7WWWnPTYzd1YUUAdNTatWuzcOHCNDQ0tK2rrKxMQ0NDFixYsNExCxYsaLd9kowZM2aT2yfJmjVr0tzc3G7pqNt/fWX6brrNpLUyeWCH1Xlh6WMdfk2AbZWgbSu9/tbr77pNRUVF3njrjS6oBoDeaPVbqwvZBoCu99JLL6WlpSV1dXXt1tfV1aWxsXGjYxobGzdr+ySZPn16amtr25b6+voO1/j6ujc28f2c9t5Y/WqHXxNgWyVo20r1NfV5z3bvecdt3mp9K/sPdKceALbMyEEj07ey7yaf71vZNyMHjezCigDobqZOnZqmpqa2Zfny5R0ee8Dg4XnrXT4Z1q5JBu8xfCurBOj9BG1bqf92/XPyyJPTp6LPRp+vSEV2qN4hJ+x3QhdXBkBv8YUPf6Ht5gcbs651XT7/4c93YUUAdNTOO++cPn36ZOXKle3Wr1y5MoMGDdromEGDBm3W9klSVVWVmpqadktHnTjh4gxYm1Rs4mo4fVqTUysOTr/+7+3wawJsqwRtBfjmkd/MsJ2HbRC29anokz6VffIf/+s/Ut23ukzVAdDTHbXHUTntw6clSbub71T+vzY+5SNTcsRuR5SjNADeRb9+/XLQQQdl3rx5betaW1szb968jB49eqNjRo8e3W77JLn77rs3uf3Wes+OA3PdPl9LZSnp8xfXauvTmuzfXJ1zp/ysU94boLcRtBWgtro2v/nMb3LWYWdlp/47JUkqKyozdu+xuf8z9+fYDxxb5goB6MkqKipy+TGX58fH/zj77rJv2/r96/bPtX93bb5z9HfKWB0A72bKlCm58sor85Of/CSPP/54Pv/5z2f16tWZPHlykmTixImZOnVq2/ann3565syZk0svvTRPPPFEzjvvvDz00EP54he/2Gk1/u2Eb+W+0Vfm2FV1qWxdv27nNyrytXwsvz776Qx435BOe2+A3qSiVCpt+naZ26jm5ubU1tamqalps065Ttbf+a3pzaZsv932qepb1UkVAvQsW3Nc7Y22dn80r2lORSoyoGpAJ1QH0PP0hD7z/e9/P5dcckkaGxszYsSIfO9738uoUaOSJEcccUR22223zJ49u237m266KWeffXaWLl2aD3zgA7n44otz7LEd/wP+1uyTN197NW+seiW1A+tT2WfT1wgF2JZ09LgqaNuIntCoAXoSx9X27A+AYjmubsg+AShWR4+rvjoKAAAAAAUQtAEAAABAAQRtAAAAAFAAQRsAAAAAFEDQBgAAAAAFELQBAAAAQAEEbQAAAABQAEEbAAAAABRA0AYAAAAABRC0AQAAAEAB+pa7gO6oVColSZqbm8tcCUDv8Pbx9O3j67ZOnwEolj6zIb0GoFgd7TWCto1YtWpVkqS+vr7MlQD0LqtWrUptbW25yyg7fQagc+gzf6bXAHSOd+s1FSV/9tlAa2trnn/++QwYMCAVFRWbPb65uTn19fVZvnx5ampqOqHCbYP9WAz7sRj249YplUpZtWpVhgwZkspKVy3QZ96dOfZ8vX1+iTl2J/rMhram1/SUf/fuzn4shv1YDPtx63W01zijbSMqKyvzV3/1V1v9OjU1NX6AC2A/FsN+LIb9uOWcYfBn+kzHmWPP19vnl5hjd6HPtFdEr+kJ/+49gf1YDPuxGPbj1ulIr/HnHgAAAAAogKANAAAAAAogaOsEVVVVmTZtWqqqqspdSo9mPxbDfiyG/Uh3si38PJpjz9fb55eYI72Xf/di2I/FsB+LYT92HTdDAAAAAIACOKMNAAAAAAogaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgbQvNnDkzu+22W6qrqzNq1Kg8+OCD77j9TTfdlGHDhqW6ujoHHHBA7rrrri6qtHvbnP04e/bsVFRUtFuqq6u7sNru5957783YsWMzZMiQVFRU5Pbbb3/XMfPnz8+HPvShVFVVZa+99srs2bM7vc7ubnP34/z58zf4WayoqEhjY2PXFMw2YVvoM5szxyuvvDIf+9jHsuOOO2bHHXdMQ0PDu+6T7mBz/x3fdv3116eioiLjxo3r3AK30ubO79VXX81pp52WwYMHp6qqKh/84Ae7/c/q5s5xxowZ2XvvvdO/f//U19fnzDPPzJtvvtlF1W4+v0tsu7aFPtMVfJ7Zeo5DxfCZpvsQtG2BG264IVOmTMm0adPy8MMPZ/jw4RkzZkxeeOGFjW5///3356STTsrJJ5+cRx55JOPGjcu4cePy+9//vosr7142dz8mSU1NTVasWNG2PPvss11YcfezevXqDB8+PDNnzuzQ9kuWLMlxxx2XI488MosWLcoZZ5yRU045JXPnzu3kSru3zd2Pb3vyySfb/TwOHDiwkypkW7Mt9JnNneP8+fNz0kkn5Ze//GUWLFiQ+vr6HH300Xnuuee6uPKO25I+lyRLly7Nl770pXzsYx/rokq3zObOb+3atTnqqKOydOnS3HzzzXnyySdz5ZVXZujQoV1cecdt7hx/+tOf5qyzzsq0adPy+OOP50c/+lFuuOGGfO1rX+viyjvO7xLbpm2hz3QFn2eK4ThUDJ9pupESm+2QQw4pnXbaaW2PW1paSkOGDClNnz59o9v//d//fem4445rt27UqFGlf/7nf+7UOru7zd2PP/7xj0u1tbVdVF3Pk6R02223veM2X/nKV0r77bdfu3Xjx48vjRkzphMr61k6sh9/+ctflpKUXnnllS6piW3PttBnNneOf2ndunWlAQMGlH7yk590VolbbUvmuG7dutKhhx5auuqqq0qTJk0qHX/88V1Q6ZbZ3PldccUVpT322KO0du3aripxq23uHE877bTSxz/+8XbrpkyZUvroRz/aqXUWxe8S245toc90BZ9niuc4VAyfacrLGW2bae3atVm4cGEaGhra1lVWVqahoSELFizY6JgFCxa02z5JxowZs8nttwVbsh+T5LXXXsuuu+6a+vr6HH/88Xnssce6otxew89isUaMGJHBgwfnqKOOym9+85tyl0MvsS30mS3tAf/T66+/nrfeeis77bRTZ5W5VbZ0jt/85jczcODAnHzyyV1R5hbbkvndcccdGT16dE477bTU1dVl//33z4UXXpiWlpauKnuzbMkcDz300CxcuLDtq2PPPPNM7rrrrhx77LFdUnNX6GnHGza0LfSZruDzTPn4eSyWzzTFE7RtppdeeiktLS2pq6trt76urm6T32VubGzcrO23BVuyH/fee+9cffXV+dnPfpbrrrsura2tOfTQQ/OnP/2pK0ruFTb1s9jc3Jw33nijTFX1PIMHD86sWbNyyy235JZbbkl9fX2OOOKIPPzww+UujV5gW+gzWzLHv/TVr341Q4YM2eAX7e5iS+Z433335Uc/+lGuvPLKrihxq2zJ/J555pncfPPNaWlpyV133ZVzzjknl156aS644IKuKHmzbckcP/3pT+eb3/xmDjvssGy33XbZc889c8QRR3Trr45uLr9L9HzbQp/pCj7PlI/jUDF8puk8fctdAHTU6NGjM3r06LbHhx56aPbZZ5/8+7//e84///wyVsa2Zu+9987ee+/d9vjQQw/N008/ne9+97u59tpry1gZbBu+/e1v5/rrr8/8+fN7zUWkV61alQkTJuTKK6/MzjvvXO5yOkVra2sGDhyYH/7wh+nTp08OOuigPPfcc7nkkksybdq0cpdXiPnz5+fCCy/MD37wg4waNSqLFy/O6aefnvPPPz/nnHNOucsDysznGboTn2k6j6BtM+28887p06dPVq5c2W79ypUrM2jQoI2OGTRo0GZtvy3Ykv34l7bbbruMHDkyixcv7owSe6VN/SzW1NSkf//+ZaqqdzjkkENy3333lbsMeoFtoc9sTQ/4zne+k29/+9u55557cuCBB3ZmmVtlc+f49NNPZ+nSpRk7dmzbutbW1iRJ37598+STT2bPPffs3KI3w5b8Gw4ePDjbbbdd+vTp07Zun332SWNjY9auXZt+/fp1as2ba0vmeM4552TChAk55ZRTkiQHHHBAVq9enc9+9rP5+te/nsrKnv9lEr9L9HzbQp/pCj7PlI/jUOfxmaYYPb/bd7F+/frloIMOyrx589rWtba2Zt68ee3+OvE/jR49ut32SXL33XdvcvttwZbsx7/U0tKSRx99NIMHD+6sMnsdP4udZ9GiRX4WKcS20Ge2tAdcfPHFOf/88zNnzpwcfPDBXVHqFtvcOQ4bNiyPPvpoFi1a1Lb87d/+bdsd1err67uy/He1Jf+GH/3oR7N48eK2ADFJ/vjHP2bw4MHdLmRLtmyOr7/++gZh2tvBYqlU6rxiu1BPO96woW2hz3QFn2fKx89j5/GZpiDlvhtDT3T99deXqqqqSrNnzy794Q9/KH32s58t7bDDDqXGxsZSqVQqTZgwoXTWWWe1bf+b3/ym1Ldv39J3vvOd0uOPP16aNm1aabvttis9+uij5ZpCt7C5+/Eb3/hGae7cuaWnn366tHDhwtKJJ55Yqq6uLj322GPlmkLZrVq1qvTII4+UHnnkkVKS0mWXXVZ65JFHSs8++2ypVCqVzjrrrNKECRPatn/mmWdK22+/fenLX/5y6fHHHy/NnDmz1KdPn9KcOXPKNYVuYXP343e/+93S7bffXnrqqadKjz76aOn0008vVVZWlu65555yTYFeZlvoM5s7x29/+9ulfv36lW6++ebSihUr2pZVq1aVawrvanPn+Je6+11HN3d+y5YtKw0YMKD0xS9+sfTkk0+Wfv7zn5cGDhxYuuCCC8o1hXe1uXOcNm1aacCAAaX/83/+T+mZZ54p/eIXvyjtueeepb//+78v1xTeld8ltk3bQp/pCj7PFMNxqBg+03QfgrYtdPnll5fe//73l/r161c65JBDSv/1X//V9tzhhx9emjRpUrvtb7zxxtIHP/jBUr9+/Ur77bdf6c477+ziirunzdmPZ5xxRtu2dXV1pWOPPbb08MMPl6Hq7uPtWzL/5fL2fps0aVLp8MMP32DMiBEjSv369SvtsccepR//+MddXnd3s7n78aKLLirtueeeperq6tJOO+1UOuKII0r/+Z//WZ7i6bW2hT6zOXPcddddN/r/02nTpnV94Zthc/8d/6fuHrSVSps/v/vvv780atSoUlVVVWmPPfYofetb3yqtW7eui6vePJszx7feeqt03nnntfWI+vr60he+8IXSK6+80vWFd5DfJbZd20Kf6Qo+z2w9x6Fi+EzTfVSUSr3kPHYAAAAAKCPXaAMAAACAAgjaAAAAAKAAgjYAAAAAKICgDQAAAAAKIGgDAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNuhhXnzxxQwaNCgXXnhh27r7778//fr1y7x588pYGQC9xTXXXJP3ve99WbNmTbv148aNy4QJE8pUFQC9hc809GYVpVKpVO4igM1z1113Zdy4cbn//vuz9957Z8SIETn++ONz2WWXlbs0AHqBN954I4MHD86VV16ZE044IUnywgsvZOjQofnFL36RI488sswVAtDT+UxDbyVogx7qtNNOyz333JODDz44jz76aH7729+mqqqq3GUB0Et84QtfyNKlS3PXXXclSS677LLMnDkzixcvTkVFRZmrA6A38JmG3kjQBj3UG2+8kf333z/Lly/PwoULc8ABB5S7JAB6kUceeSQf/vCH8+yzz2bo0KE58MADc8IJJ+Scc84pd2kA9BI+09AbuUYb9FBPP/10nn/++bS2tmbp0qXlLgeAXmbkyJEZPnx4rrnmmixcuDCPPfZY/umf/qncZQHQi/hMQ2/kjDbogdauXZtDDjkkI0aMyN57750ZM2bk0UcfzcCBA8tdGgC9yBVXXJEZM2bkqKOOylNPPZW5c+eWuyQAegmfaeitBG3QA335y1/OzTffnN/97nd573vfm8MPPzy1tbX5+c9/Xu7SAOhFmpqaMmTIkKxbty7XXHNNxo8fX+6SAOglfKaht/LVUehh5s+fnxkzZuTaa69NTU1NKisrc+211+bXv/51rrjiinKXB0AvUltbm09+8pN573vfm3HjxpW7HAB6CZ9p6M2c0QYAwCZ94hOfyH777Zfvfe975S4FAKDbE7QBALCBV155JfPnz8+nPvWp/OEPf8jee+9d7pIAALq9vuUuAACA7mfkyJF55ZVXctFFFwnZAAA6yBltAAAAAFAAN0MAAAAAgAII2gAAAACgAII2AAAAACiAoA0AAAAACiBoAwAAAIACCNoAAAAAoACCNgAAAAAogKANAAAAAAogaAMAAACAAvz/e/gHozxT/KgAAAAASUVORK5CYII=", + "text/plain": [ + "
    " + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], "source": [ "plt.figure(figsize=(15,5))\n", "plt.subplot(131)\n", @@ -615,9 +697,18 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 13, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "-0.06089330922309355\n", + "-0.06055651225519789\n" + ] + } + ], "source": [ "print(calculate_exchange_tensor(pairs[0])[0]) # isotropic should be -82 meV\n", "print(calculate_exchange_tensor(pairs[1])[0]) # these should all be around -41.9 in the isotropic part\n", @@ -628,9 +719,18 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 14, "metadata": {}, - "outputs": [], + "outputs": [ + { + "ename": "SyntaxError", + "evalue": "invalid syntax (659628047.py, line 1)", + "output_type": "error", + "traceback": [ + "\u001b[0;36m Cell \u001b[0;32mIn[14], line 1\u001b[0;36m\u001b[0m\n\u001b[0;31m These are reasonably converged:\u001b[0m\n\u001b[0m ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m invalid syntax\n" + ] + } + ], "source": [ "These are reasonably converged:\n", "\n", @@ -647,10 +747,6 @@ "metadata": {}, "outputs": [], "source": [ - "# symmetrizing Hamiltonian and overlap matrix to make them hermitian \n", - "# Check if exchange field has scalar part\n", - "# parallel over integrals\n", - "\n", "========================================\n", " \n", "Atom Angstrom\n", diff --git a/test.py b/test.py deleted file mode 100644 index 766b8cd..0000000 --- a/test.py +++ /dev/null @@ -1,515 +0,0 @@ -import os -from sys import stdout -from timeit import default_timer as timer - -from tqdm import tqdm - -os.environ["OMP_NUM_THREADS"] = "1" # export OMP_NUM_THREADS=4 -os.environ["OPENBLAS_NUM_THREADS"] = "1" # export OPENBLAS_NUM_THREADS=4 -os.environ["MKL_NUM_THREADS"] = "1" # export MKL_NUM_THREADS=6 -os.environ["VECLIB_MAXIMUM_THREADS"] = "1" # export VECLIB_MAXIMUM_THREADS=4 -os.environ["NUMEXPR_NUM_THREADS"] = "1" # export NUMEXPR_NUM_THREADS=6 - -import warnings - -import numpy as np -import sisl -from mpi4py import MPI -from numpy.linalg import inv - -from grogu.useful import * - -start_time = timer() - -# this cell mimicks an input file -fdf = sisl.get_sile("./lat3_791/Fe3GeTe2.fdf") -# this information needs to be given at the input!! -scf_xcf_orientation = np.array([0, 0, 1]) # z -# list of reference directions for around which we calculate the derivatives -# o is the quantization axis, v and w are two axes perpendicular to it -# at this moment the user has to supply o,v,w on the input. -# we can have some default for this -ref_xcf_orientations = [ - dict(o=np.array([1, 0, 0]), vw=[np.array([0, 1, 0]), np.array([0, 0, 1])]), - dict(o=np.array([0, 1, 0]), vw=[np.array([1, 0, 0]), np.array([0, 0, 1])]), - dict(o=np.array([0, 0, 1]), vw=[np.array([1, 0, 0]), np.array([0, 1, 0])]), -] - -# human readable definition of magnetic entities -# magnetic_entities = [ -# dict(atom=0, ), -# dict(atom=1, ), -# dict(atom=2, ), -# dict(atom=3, l=2), -# dict(atom=4, l=2), -# dict(atom=5, l=2), -# ] -# pairs = [ -# dict(ai=3, aj=4, Ruc=np.array([0, 0, 0])), # isotropic should be -82 meV -# dict(ai=3, aj=5, Ruc=np.array([0, 0, 0])), # these should all be around -41.9 in the isotropic part -# dict(ai=4, aj=5, Ruc=np.array([0, 0, 0])), -# dict(ai=3, aj=0, Ruc=np.array([0, 0, 0])), -# dict(ai=3, aj=1, Ruc=np.array([0, 0, 0])), -# dict(ai=3, aj=2, Ruc=np.array([0, 0, 0])), -# ] -magnetic_entities = [ - dict(atom=3, l=2), - dict(atom=4, l=2), - dict(atom=5, l=2), -] - -# pair information -pairs = [ - dict(ai=0, aj=1, Ruc=np.array([0, 0, 0])), # isotropic should be -82 meV - dict( - ai=0, aj=2, Ruc=np.array([0, 0, 0]) - ), # these should all be around -41.9 in the isotropic part - dict(ai=1, aj=2, Ruc=np.array([0, 0, 0])), - dict(ai=0, aj=1, Ruc=np.array([-1, 0, 0])), - dict(ai=0, aj=2, Ruc=np.array([-1, 0, 0])), - dict(ai=0, aj=1, Ruc=np.array([1, 0, 0])), - dict(ai=0, aj=2, Ruc=np.array([1, 0, 0])), - dict(ai=0, aj=1, Ruc=np.array([0, -1, 0])), - dict(ai=0, aj=2, Ruc=np.array([0, -1, 0])), - dict(ai=0, aj=1, Ruc=np.array([0, 1, 0])), - dict(ai=0, aj=2, Ruc=np.array([0, 1, 0])), - dict(ai=1, aj=2, Ruc=np.array([-1, 0, 0])), -] - -# Brilloun zone sampling and Green function contour integral -kset = 20 -kdirs = "xy" -ebot = -30 -eset = 50 -esetp = 1000 - - -# MPI parameters -comm = MPI.COMM_WORLD -size = comm.Get_size() -rank = comm.Get_rank() -root_node = 0 -if rank == root_node: - print("Number of nodes in the parallel cluster: ", size) - -simulation_parameters = dict( - path="Not yet specified.", - scf_xcf_orientation=scf_xcf_orientation, - ref_xcf_orientations=ref_xcf_orientations, - kset=kset, - kdirs=kdirs, - ebot=ebot, - eset=eset, - esetp=esetp, - parallel_size=size, -) - -# digestion of the input -# read in hamiltonian -dh = fdf.read_hamiltonian() -try: - simulation_parameters["geom"] = fdf.read_geometry() -except: - print("Error reading geometry.") - -# unit cell index -uc_in_sc_idx = dh.lattice.sc_index([0, 0, 0]) - -setup_time = timer() - -NO = dh.no # shorthand for number of orbitals in the unit cell - -# preprocessing Hamiltonian and overlap matrix elements -h11 = dh.tocsr(dh.M11r) -h11 += dh.tocsr(dh.M11i) * 1.0j -h11 = h11.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") - -h22 = dh.tocsr(dh.M22r) -h22 += dh.tocsr(dh.M22i) * 1.0j -h22 = h22.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") - -h12 = dh.tocsr(dh.M12r) -h12 += dh.tocsr(dh.M12i) * 1.0j -h12 = h12.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") - -h21 = dh.tocsr(dh.M21r) -h21 += dh.tocsr(dh.M21i) * 1.0j -h21 = h21.toarray().reshape(NO, dh.n_s, NO).transpose(0, 2, 1).astype("complex128") - -sov = ( - dh.tocsr(dh.S_idx) - .toarray() - .reshape(NO, dh.n_s, NO) - .transpose(0, 2, 1) - .astype("complex128") -) - - -# Reorganization of Hamiltonian and overlap matrix elements to SPIN BOX representation -U = np.vstack( - [np.kron(np.eye(NO, dtype=int), [1, 0]), np.kron(np.eye(NO, dtype=int), [0, 1])] -) -# This is the permutation that transforms ud1ud2 to u12d12 -# That is this transforms FROM SPIN BOX to ORBITAL BOX => U -# the inverse transformation is U.T u12d12 to ud1ud2 -# That is FROM ORBITAL BOX to SPIN BOX => U.T - -# From now on everything is in SPIN BOX!! -hh, ss = np.array( - [ - U.T @ np.block([[h11[:, :, i], h12[:, :, i]], [h21[:, :, i], h22[:, :, i]]]) @ U - for i in range(dh.lattice.nsc.prod()) - ] -), np.array( - [ - U.T - @ np.block([[sov[:, :, i], sov[:, :, i] * 0], [sov[:, :, i] * 0, sov[:, :, i]]]) - @ U - for i in range(dh.lattice.nsc.prod()) - ] -) - - -# symmetrizing Hamiltonian and overlap matrix to make them hermitian -for i in range(dh.lattice.sc_off.shape[0]): - j = dh.lattice.sc_index(-dh.lattice.sc_off[i]) - h1, h1d = hh[i], hh[j] - hh[i], hh[j] = (h1 + h1d.T.conj()) / 2, (h1d + h1.T.conj()) / 2 - s1, s1d = ss[i], ss[j] - ss[i], ss[j] = (s1 + s1d.T.conj()) / 2, (s1d + s1.T.conj()) / 2 - -# identifying TRS and TRB parts of the Hamiltonian -TAUY = np.kron(np.eye(NO), tau_y) -hTR = np.array([TAUY @ hh[i].conj() @ TAUY for i in range(dh.lattice.nsc.prod())]) -hTRS = (hh + hTR) / 2 -hTRB = (hh - hTR) / 2 - -# extracting the exchange field -traced = [spin_tracer(hTRB[i]) for i in range(dh.lattice.nsc.prod())] # equation 77 -XCF = np.array( - [ - np.array([f["x"] for f in traced]), - np.array([f["y"] for f in traced]), - np.array([f["z"] for f in traced]), - ] -) # equation 77 - -# Check if exchange field has scalar part -max_xcfs = abs(np.array(np.array([f["c"] for f in traced]))).max() -if max_xcfs > 1e-12: - warnings.warn( - f"Exchange field has non negligible scalar part. Largest value is {max_xcfs}" - ) - -H_and_XCF_time = timer() - -# for every site we have to store 3 Greens function (and the associated _tmp-s) in the 3 reference directions -for i, mag_ent in enumerate(magnetic_entities): - parsed = parse_magnetic_entity(dh, **mag_ent) # parse orbital indexes - magnetic_entities[i]["orbital_indeces"] = parsed - magnetic_entities[i]["spin_box_indeces"] = blow_up_orbindx( - parsed - ) # calculate spin box indexes - spin_box_shape = len( - mag_ent["spin_box_indeces"] - ) # calculate size for Greens function generation - - mag_ent["energies"] = [] # we will store the second order energy derivations here - - mag_ent["Gii"] = [] # Greens function - mag_ent["Gii_tmp"] = [] # Greens function for parallelization - mag_ent["Vu1"] = [ - list([]) for _ in range(len(ref_xcf_orientations)) - ] # These will be the perturbed potentials from eq. 100 - mag_ent["Vu2"] = [list([]) for _ in range(len(ref_xcf_orientations))] - for i in ref_xcf_orientations: - mag_ent["Gii"].append( - np.zeros((eset, spin_box_shape, spin_box_shape), dtype="complex128") - ) # Greens functions for every quantization axis - mag_ent["Gii_tmp"].append( - np.zeros((eset, spin_box_shape, spin_box_shape), dtype="complex128") - ) - -# for every site we have to store 2x3 Greens function (and the associated _tmp-s) -# in the 3 reference directions, because G_ij and G_ji are both needed -for pair in pairs: - spin_box_shape_i, spin_box_shape_j = len( - magnetic_entities[pair["ai"]]["spin_box_indeces"] - ), len( - magnetic_entities[pair["aj"]]["spin_box_indeces"] - ) # calculate size for Greens function generation - - pair["energies"] = [] # we will store the second order energy derivations here - - pair["Gij"] = [] # Greens function - pair["Gji"] = [] - pair["Gij_tmp"] = [] # Greens function for parallelization - pair["Gji_tmp"] = [] - - pair["Vij"] = [ - list([]) for _ in range(len(ref_xcf_orientations)) - ] # These will be the perturbed potentials from eq. 100 - pair["Vji"] = [list([]) for _ in range(len(ref_xcf_orientations))] - - for i in ref_xcf_orientations: - pair["Gij"].append( - np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype="complex128") - ) - pair["Gij_tmp"].append( - np.zeros((eset, spin_box_shape_i, spin_box_shape_j), dtype="complex128") - ) # Greens functions for every quantization axis - pair["Gji"].append( - np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype="complex128") - ) - pair["Gji_tmp"].append( - np.zeros((eset, spin_box_shape_j, spin_box_shape_i), dtype="complex128") - ) - -site_and_pair_dictionaries_time = timer() - -kset = make_kset(dirs=kdirs, NUMK=kset) # generate k space sampling -wkset = np.ones(len(kset)) / len(kset) # generate weights for k points -kpcs = np.array_split(kset, size) # split the k points based on MPI size -kpcs[root_node] = tqdm(kpcs[root_node], desc="k loop", file=stdout) - -k_set_time = timer() - -# this will contain all the data needed to calculate the energy variations upon rotation -hamiltonians = [] - -# iterate over the reference directions (quantization axes) -for i, orient in enumerate(ref_xcf_orientations): - # obtain rotated exchange field - R = RotMa2b(scf_xcf_orientation, orient["o"]) - rot_XCF = np.einsum("ij,jklm->iklm", R, XCF) - rot_H_XCF = sum( - [np.kron(rot_XCF[i], tau) for i, tau in enumerate([tau_x, tau_y, tau_z])] - ) - rot_H_XCF_uc = rot_H_XCF[uc_in_sc_idx] - - # obtain total Hamiltonian with the rotated exchange field - rot_H = hTRS + rot_H_XCF # equation 76 - - hamiltonians.append( - dict(orient=orient["o"], H=rot_H, rotations=[]) - ) # store orientation and rotated Hamiltonian - - for u in orient[ - "vw" - ]: # these are the infinitezimal rotations (for now) perpendicular to the quantization axis - Tu = np.kron(np.eye(NO, dtype=int), tau_u(u)) # section 2.H - - Vu1 = 1j / 2 * commutator(rot_H_XCF_uc, Tu) # equation 100 - Vu2 = 1 / 8 * commutator(commutator(Tu, rot_H_XCF_uc), Tu) # equation 100 - - for mag_ent in magnetic_entities: - mag_ent["Vu1"][i].append( - Vu1[:, mag_ent["spin_box_indeces"]][mag_ent["spin_box_indeces"], :] - ) # fill up the perturbed potentials (for now) based on the on-site projections - mag_ent["Vu2"][i].append( - Vu2[:, mag_ent["spin_box_indeces"]][mag_ent["spin_box_indeces"], :] - ) - - for pair in pairs: - ai = magnetic_entities[pair["ai"]][ - "spin_box_indeces" - ] # get the pair orbital sizes from the magnetic entities - aj = magnetic_entities[pair["aj"]]["spin_box_indeces"] - pair["Vij"][i].append( - Vu1[:, ai][aj, :] - ) # fill up the perturbed potentials (for now) based on the on-site projections - pair["Vji"][i].append(Vu1[:, aj][ai, :]) - -reference_rotations_time = timer() - -if rank == root_node: - print("Number of magnetic entities being calculated: ", len(magnetic_entities)) - print( - "We have to calculate the Greens function for three reference direction and we are going to calculate 15 energy integrals per site." - ) - print(f"The shape of the Hamiltonian and the Greens function is {NO}x{NO}.") -comm.Barrier() -# ---------------------------------------------------------------------- - -# make energy contour -# we are working in eV now ! -# and sisil shifts E_F to 0 ! -cont = make_contour(emin=ebot, enum=eset, p=esetp) -eran = cont.ze - -# ---------------------------------------------------------------------- -# sampling the integrand on the contour and the BZ -for k in kpcs[rank]: - wk = wkset[rank] # weight of k point in BZ integral - for i, hamiltonian_orientation in enumerate( - hamiltonians - ): # iterate over reference directions - # calculate Greens function - H = hamiltonian_orientation["H"] - HK, SK = hsk(H, ss, dh.sc_off, k) - Gk = inv(SK * eran.reshape(eset, 1, 1) - HK) - - # store the Greens function slice of the magnetic entities (for now) based on the on-site projections - for mag_ent in magnetic_entities: - mag_ent["Gii_tmp"][i] += ( - Gk[:, mag_ent["spin_box_indeces"]][..., mag_ent["spin_box_indeces"]] - * wk - ) - - for pair in pairs: - # add phase shift based on the cell difference - phase = np.exp(1j * 2 * np.pi * k @ pair["Ruc"].T) - - # get the pair orbital sizes from the magnetic entities - ai = magnetic_entities[pair["ai"]]["spin_box_indeces"] - aj = magnetic_entities[pair["aj"]]["spin_box_indeces"] - - # store the Greens function slice of the magnetic entities (for now) based on the on-site projections - pair["Gij_tmp"][i] += Gk[:, ai][..., aj] * phase * wk - pair["Gji_tmp"][i] += Gk[:, aj][..., ai] * phase * wk - -# summ reduce partial results of mpi nodes -for i in range(len(hamiltonians)): - for mag_ent in magnetic_entities: - comm.Reduce(mag_ent["Gii_tmp"][i], mag_ent["Gii"][i], root=root_node) - - for pair in pairs: - comm.Reduce(pair["Gij_tmp"][i], pair["Gij"][i], root=root_node) - comm.Reduce(pair["Gji_tmp"][i], pair["Gji"][i], root=root_node) - -green_function_inversion_time = timer() - -if rank == root_node: - # iterate over the magnetic entities - for tracker, mag_ent in enumerate(magnetic_entities): - # iterate over the quantization axes - for i, Gii in enumerate(mag_ent["Gii"]): - storage = [] - # iterate over the first and second order local perturbations - for Vu1, Vu2 in zip(mag_ent["Vu1"][i], mag_ent["Vu2"][i]): - # The Szunyogh-Lichtenstein formula - traced = np.trace((Vu2 @ Gii + 0.5 * Gii @ Vu1 @ Gii), axis1=1, axis2=2) - # evaluation of the contour integral - storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we))) - - # fill up the magnetic entities dictionary with the energies - mag_ent["energies"].append(storage) - - # iterate over the pairs - for tracker, pair in enumerate(pairs): - # iterate over the quantization axes - for i, (Gij, Gji) in enumerate(zip(pair["Gij"], pair["Gji"])): - site_i = magnetic_entities[pair["ai"]] - site_j = magnetic_entities[pair["aj"]] - - storage = [] - # iterate over the first order local perturbations in all possible orientations for the two sites - for Vui in site_i["Vu1"][i]: - for Vuj in site_j["Vu1"][i]: - # The Szunyogh-Lichtenstein formula - traced = np.trace((Vui @ Gij @ Vuj @ Gji), axis1=1, axis2=2) - # evaluation of the contour integral - storage.append(np.trapz(-1 / np.pi * np.imag(traced * cont.we))) - - # fill up the pairs dictionary with the energies - pairs[tracker]["energies"].append(storage) - - end_time = timer() - - print( - "############################### GROGU OUTPUT ###################################" - ) - print( - "================================================================================" - ) - print("Input file: ") - print(simulation_parameters["path"]) - print( - "Number of nodes in the parallel cluster: ", - simulation_parameters["parallel_size"], - ) - print( - "================================================================================" - ) - try: - print("Cell [Ang]: ") - print(simulation_parameters["geom"].cell) - except: - print("Geometry could not be read.") - print( - "================================================================================" - ) - print("DFT axis: ") - print(simulation_parameters["scf_xcf_orientation"]) - print("Quantization axis and perpendicular rotation directions:") - for ref in ref_xcf_orientations: - print(ref["o"], " --» ", ref["vw"]) - print( - "================================================================================" - ) - print("number of k points: ", simulation_parameters["kset"]) - print("k point directions: ", simulation_parameters["kdirs"]) - print( - "================================================================================" - ) - print("Parameters for the contour integral:") - print("Ebot: ", simulation_parameters["ebot"]) - print("Eset: ", simulation_parameters["eset"]) - print("Esetp: ", simulation_parameters["esetp"]) - print( - "================================================================================" - ) - print("Atomic informations: ") - print("") - print("") - print("Not yet specified.") - print("") - print("") - print( - "================================================================================" - ) - print("Exchange [meV]") - print( - "--------------------------------------------------------------------------------" - ) - print("Atom1 Atom2 [i j k] d [Ang]") - print( - "--------------------------------------------------------------------------------" - ) - for pair in pairs: - J_iso, J_S, D = calculate_exchange_tensor(pair) - J_iso = J_iso * sisl.unit_convert("eV", "meV") - J_S = J_S * sisl.unit_convert("eV", "meV") - D = D * sisl.unit_convert("eV", "meV") - - print(print_atomic_indices(pair, magnetic_entities, dh)) - print("Isotropic: ", J_iso) - print("DMI: ", D) - print("Symmetric-anisotropy: ", J_S) - print("") - - print( - "================================================================================" - ) - print("Runtime information: ") - print("Total runtime: ", end_time - start_time) - print( - "--------------------------------------------------------------------------------" - ) - print("Initial setup: ", setup_time - start_time) - print( - f"Hamiltonian conversion and XC field extraction: {H_and_XCF_time - setup_time:.3f} s" - ) - print( - f"Pair and site datastructure creatrions: {site_and_pair_dictionaries_time - H_and_XCF_time:.3f} s" - ) - print( - f"k set cration and distribution: {k_set_time - site_and_pair_dictionaries_time:.3f} s" - ) - print(f"Rotating XC potential: {reference_rotations_time - k_set_time:.3f} s") - print( - f"Greens function inversion: {green_function_inversion_time - reference_rotations_time:.3f} s" - ) - print( - f"Calculate energies and magnetic components: {end_time - green_function_inversion_time:.3f} s" - )