From 3de31e212333fd154abeb89602b6a437cdba773c Mon Sep 17 00:00:00 2001 From: NehZio Date: Thu, 3 Mar 2022 16:38:58 +0100 Subject: [PATCH] Revert "Added plane orientation" This reverts commit 66a287108621ff549a0441345a88575022033439. --- crystal_mec.py | 42 +-- src/__pycache__/out.cpython-38.pyc | Bin 5308 -> 0 bytes src/__pycache__/out.cpython-39.pyc | Bin 5314 -> 0 bytes src/__pycache__/read_input.cpython-38.pyc | Bin 5107 -> 0 bytes src/__pycache__/read_input.cpython-39.pyc | Bin 5068 -> 0 bytes src/__pycache__/utils.cpython-38.pyc | Bin 8091 -> 0 bytes src/__pycache__/utils.cpython-39.pyc | Bin 8202 -> 0 bytes src/read_input.py | 26 +- src/read_input.py~ | 196 ----------- src/utils.py | 10 +- src/utils.py~ | 390 ---------------------- 11 files changed, 6 insertions(+), 658 deletions(-) delete mode 100644 src/__pycache__/out.cpython-38.pyc delete mode 100644 src/__pycache__/out.cpython-39.pyc delete mode 100644 src/__pycache__/read_input.cpython-38.pyc delete mode 100644 src/__pycache__/read_input.cpython-39.pyc delete mode 100644 src/__pycache__/utils.cpython-38.pyc delete mode 100644 src/__pycache__/utils.cpython-39.pyc delete mode 100644 src/read_input.py~ delete mode 100644 src/utils.py~ diff --git a/crystal_mec.py b/crystal_mec.py index e4fcadb..42c7fc8 100644 --- a/crystal_mec.py +++ b/crystal_mec.py @@ -31,7 +31,7 @@ if __name__=='__main__': inputFile = sys.argv[1] # Reads all the parameters from the input file - rB , rPP, center, X, Y, Z, xOy, yOz, xOz, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translation = read_input(inputFile) + rB , rPP, center, X, Y, Z, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translation = read_input(inputFile) if verbose > 0: out_input_param(rB , rPP, center, X, Y, Z, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translation) @@ -79,59 +79,25 @@ if __name__=='__main__': # Orienting the big cell - if xOy != []: - a = find_center(xOy[0], coordinates) - b = a - w = [a] - while np.absolute ( np.absolute( np.dot( a / np.linalg.norm(a) , b / np.linalg.norm(a) ) ) - 1 ) <= 1e-6: - b = find_center(xOy[1], coordinates, without=w) - w.append(b) - c = np.cross(a, b) - k = [0,0,1] - M = rotation_matrix(c, k) - coordinates = rotate(M, coordinates) - if xOz != []: - a = find_center(xOz[0], coordinates) - b = a - w = [a] - while np.absolute ( np.absolute( np.dot( a / np.linalg.norm(a) , b / np.linalg.norm(a) ) ) - 1 ) <= 1e-6: - b = find_center(xOz[1], coordinates, without=w) - w.append(b) - c = np.cross(a, b) - k = [0,1,0] - M = rotation_matrix(c, k) - coordinates = rotate(M, coordinates) - if yOz != []: - a = find_center(yOz[0], coordinates) - b = a - w = [a] - while np.absolute ( np.absolute( np.dot( a / np.linalg.norm(a) , b / np.linalg.norm(a) ) ) - 1 ) <= 1e-6: - b = find_center(yOz[1], coordinates, without=w) - w.append(b) - c = np.cross(a, b) - k = [1,0,0] - M = rotation_matrix(c, k) - coordinates = rotate(M, coordinates) if X != []: k = [1,0,0] xVec = find_center(X,coordinates) - M = rotation_matrix(xVec, k) + M = rotation_matrix(k,xVec) coordinates = rotate(M, coordinates) - if Y != []: k = [0,1,0] yVec = find_center(Y,coordinates) - M = rotation_matrix(yVec, k) + M = rotation_matrix(k,yVec) coordinates = rotate(M, coordinates) if Z != []: k = [0,0,1] zVec = find_center(Z,coordinates) - M = rotation_matrix(zVec, k) + M = rotation_matrix(k,zVec) coordinates = rotate(M, coordinates) diff --git a/src/__pycache__/out.cpython-38.pyc b/src/__pycache__/out.cpython-38.pyc deleted file mode 100644 index 0caae17558dfc52a3b08c27cf76c9b49a3b602a7..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 5308 zcma)AO>Er873T0~f2beZiY)8LN!Zj%wTiuxoJet1$BiWaG^s7SjZ!$<#s;-Rtr&61 z)sS0R1Q#f(0L>+U3-nav;N+H5b16{t)KgLP+!N7TiUK~ir=X{PZ^)%3xo*M*hi~4z z_vX#aH*enTyF){!g6HqQ`;)CqD9V>qIs9m-yoonlMZp!$Jf*?pS8b^Hs-D)?8#+_C z#`O=ChQWt;nO8n!jS{bN6Pz+1;Dh*9ID1E_4L^n*WxdAI3U|ZE_M1*~l34k37T*fq zHN4@Bt|$@uOv+=*0b@j~iBiy7!U|f?+l^$f;FXf0f-Wb+1zkx_6m&H?Sy zqY+Cg2aFHzsG|x$UDWhZ<&2Wo#`u||{C>}on(+@9`t@bM<4HA{oKW}#Ka*3#JJQyc znvYH8L!?Dg==a2d2HkpZ3f@d++++$Pl()&aF!h40kC}+)t2eZ<9>+5C>J^vHltFSw?qd!!|+wkLC(uxWzW75h? z$>bq_)!p-va-aY`EYS+RUqbKOQjf-;>76^K#|u56EAQ`*(3_Te`ZM)DDfEm&@88gS zS?U?j)caeemrU^Q(66QKXb&}EL30K?Lf?c$-vmE|K8&iA_o9irkgdsHN<}ZTM@FUX zF6X^G2iseae?`iVq86cczNpa*U7#9zI9>GcS8(h3ni4lhSB9lJ-Mw z7f8?Ne=N77i*v9upW7ksQ9EO>b57b>cp4FX{J`>8nPsjZvRqt`U%cIoI^D?Hay`d# zL+jA9Y)yw|JiTg1+m^7o+da|>K?qn2OZA2MEr{I-oh}b9bb`q7BiHtf{!A~-OBWA+1w*KPEDIG>eAZlv zg0|bVIzi}W+SyQ0F2^&Hp1zowk>koZ$H~WZa>S$N#pDn!e~-)0f~&s=8rA zvxSqEJD5R7_$dqH8YHc)!1IDV*Kc7?WeaBB^DMj5acq%USPWA=a@$Tkfbqa2iqK3I zoG6;(^LXTVkBebLw=tsljhptSe3o{#&_1U8-ZxM_+Fj%8HCk9<`NV5S z*DeyaA9@mc+4^ldw1l4PaJkPG>#t(=Xcp<46Y=>C;Xo%J?h_>ZZ#>JP0#%k3|CD;~ zh@gdCnDr+dXSeHm?k4c>$k}SjB@Zn-l2y-byHQB%PE#3n2&AFa4XO8QKmPy=ge|lk za5H}EDdQb1WM*Yc*sV5})@lZU;I5C(Ld(q;mts&~p#ES#Hsibh$J8bTre~N5Eovnj zmJK{o<$g2pfK`WQmT!>dCo?v^NI1d1@q-A6gIzy#Iflsl-gRVvdlKKmZ3p7ZaXCj{ zyz~rwX{EpPsTs(~ejSr4%D`M&@%y2}pv{ArKv1k>%KfsN-!3FMEi>k%3me__G=Ypu?l`KvsHGIonzx{9u#qA*s%PL zveMU3C)o>!Budcgzr*Y-YA2-yqF!R+JbXzt+w)S@51IzKmOn%IR`F8I4^5mUk&-aZ z4%D9UD`ksum1`eX66Go~f?g@nxxS-{S#I>ozRpX%a>PDRlM>w@C)!R){P`R9LUQF#72Js^ zFFCQEvb`pqulaPAy1}2H9nwqNhc`*aQ&sl1C zbc)Fa@d{X}dD#QBL|V9h)wm&UeL6?@Ct~&j;KRk&|*F& z%?+HRL&Ax=!iQ!TpcNOuj*V&8GR=54Kl!nx7v{5VE(7?}6;tLFf)F9daIKQ+NJxBO zr<0btzPsCXXwx>MREynCYDjMdAt^PQUJyF9p>*Jo^Nk9nb?Zp9#6<`-jJ<8QiLFeL zP^4@tWnM$W< z5jG5J1RG&i)lm-1Z8OMbR6$qOCJnLw(2aTkFExaJK!M$WTP${iPCJ3?-1tQ689;t1 zq1(W3dgVla#CSO|2){Z0h*OamC$#cjCBdx?H@hm(SjKKx#f~w7+)~aTm{20@q62mY zXqNT}-YA6bHM#DwM(Oq~t4O(klI(fVCk}u24_9(JE_hkG zW{TI~N(Lx`YrHVMlxLw?-cIuU9R3)cCE%t6V#>^#CYH$7_fVu-7>SxG#)$C-G0J(a zkXBuvI}dJRr-<)U1LI;vej|w{@g^!6oG9>Ta2m|ulpR+-*x@(|AW0nrihv+sio$Az zjj|;sZldi-Y@jGw{Qpr$)MzO^7SSCHF^{uxlt9B@dMYrba;$aMazMA=FM3)6?u63Q z@zwA(fCxP?fV>h{fwB_Uc68C;IwCJm$kILF3U~w-DjX;Qu6~&qd7V)`{H$n~b#bp~ z&*bk#+9N$FK_A#1N1WxA9Yx%uK6dXVWigwuM@DWVvs`Lv5@WNwqlo1!DnEiQ^r|1L z7%ToYlDEz(&;ElUFVA7TRs08|l8-cQS)#w1zC}`f)ryR$0N}fbNVW< zL4RQZtNp+zTaH%f5&(n^!vh)L1RVh4brd*FH^a0<*BhUwr7#eY!_!K3aYIxy5y0Cm zhi*b*5c0LEglJ>1@1q1KoeZKz1^tsIrl@!uMWbXln@$H=(av$0PGreOG3fh!gW&-b zv*^ckb(D<(sTF3*m{$URe?rlaomD9=(ofKtQ&Z(3A&-(;H8uQhyR$C`C@7ewRiruf vZsdAlHgrlPWZ0DWmnIH5_+&FxF1G{T^_(kodWJLs=n@;L@`#3ffl zZe$T$pojw0mjEu%Q;~y{TP|%b1&W?}D)I;PMD&)Tg^z84qUfo=H{{ZiTsPr@GjAT> zd-G=IH*emoH85Z(c>esGKUm7RqI^M(!;gl>ckl))D7eCztJImS)w+sWb+wjW*O|gK zu79A^3w(f=c=-cXFY*dEz$x)Q-jBM>*&9lA@CocF>s6MNxf6tz*RUHC#LAv2)MdP@ zc!Tr0qJ-!(E{!M$j1jHIN=|Dr%V|CPUWog1UNIiX=~6tH)8%+5rz`P^oHpW*D^i1nBXuZIP zLl&107$4kGhZTN0Z|TFzsFJlt_~=o7ujf!r`3DUBda~cKxDro{D}0=fX4K%0w6&#X zV^jG6X^|BA9XimUTkB53n~9VgPhy1fdpyYpM*$1c<&%8sLzO9y*e=IA1Nj%F_r?+L z(^k)WYMqhZAL8fUEBqxnYX&@+mEN0QHf!klukc=l-GLqbkt)6iKfWfdsIW33t(=V~ z4*9F@UI>)~1?XY1mg~I_y{}6>8h@&H?wB6W^?!G6`q0Fikgdu-ig_QiN5-V> zUe5Y>0k$zd{<4%GMoXY|A#V|;FH#FVoX&gr9Nb#AqQuRxl_D!Xjb2a56_}%5N&Atu z`|Er*Mlw6PI14-TnH}ODwKD=c=cJv5XYtTO4=n#HwagX7l}qc4NjXv$a) z{g%@(+kW7r+UZbGUWsNTJ$)%PBgd6;j+0O5+=zyajcr@l=JYFpIep17r!T{oM0J9? zW(hkfwK0RX@Ddh8RY;m!zU%sXj@QJT$}gCC*EOwn+qOh%VKGSb&}rFGAI1ZdC;}r< zaGq$k$D>ondt40ax`h!%uidgXZI`@V-L^#2CThLY+O&n~Z<(I$G`BbXj$Fgr)&r;2 zX_=e8$4%Gvn&Eczs_coBu%3&xOIwo63+YPUN{MN6srF6ask#an*@A|Jb=9Zp z%rv8)Q15NwH?bqr{)BDsb{yB)#9pxF%r)fl2c{Lurt7qvFdz)jv<7VgY+!Z*>iznA z?_#mAowj{$M6W+%yp2Uot!xRa*}^iL4c`~s@z7adI@v-~Jn9M59_&X(^w$3|wMjwh zIc5TrzLE{o0ye2}zu~*Ut<6)*H_7tT8Jk`t%wP|BehBozt`|5QL*za0+A`QZjdj7c z4e={cDMMnk^c;(@WeS5g2- zns#+6QN>E4iMw}`Vnc43guRupw-Z)R%Dtse3|~h18<8q*9#*B+DM2GlvpOr zO4w$?8c9LM(?q8jniOQ%gt={hk4}T6VAC3>j2m!1>BCID?xjaSQjrbG?}rCBZI5EA zFOo_#C;JXzJYi14-u+Y&FT=C{EX{5EEql(j{iYS*l;Cq~expM$3FcOPZqKcXeWbwd z#oO0c=K|4~!_aE&{bUIBeaFLDc3++?MHmYnqEJ}B{FSr1*D5oUmGL{zR@e!4j*YQ- zP{f&GgR&lG#V@7K9FizPtM?AFDYQ;V3q+k|;v#%WG|P1p)$9R z^s$255#=aD>j~Ry&{^Ai=*Yz#M|lBnFpQ#!alF&tWig|^exdH_Z=_(7h^nfITALSD z2!6_vf=efw{2(rXl^9oDAiv?a+SfMF{wH1~n%)ymXd_DGajACI^&6HOTm#XQ7juvm zw2)88pADR+1HukF!h>c9fE6_|R+x57!-!_H3;^4CVLsjBGI&4RFl2rq2n~V^SIdcx z#KZ$`+DWnFIlCR3Hf$qIw8&{E1?jCAglJM|xPD+)2a>)+&ezM7-mN3y60b<1y=|v~ zZA?*6By1~TZe7E~h%2cHa=u=cGqvV4Vmc@&Cjuv03w zgQ}~i*dVA=*aWkxj&e|LnSM5-iW`t@&=7kM-K+cX%A4%(Kwu}}E{mO@(=OmTFMOnR z3qXA_raQqeyQNrv%y=m-5Oy={5vMHiO-SY4a*Vqh?spX+v4owlg56>OwZ)9zH?D-( zK?iIUNS5{p+o?_wIw@u36zvR^EN&IDov|Xd2IVd}li^m(!~c5u#t)Zo-oC!RQ6IYd z#`Pbq-MM>b1FhGWZ`Jjc<&Ezmtg}0Hb=6JSs$BO-qqO_BnWtVrN%}nK5r^yj!;OrN za$cIY8Da^pq<|v0MhnwRSvH#H?Slo;P8 zMk&h^l8WPT`@t>j6hZ%iB#etSSwn%OX?;^TQP5A})StpBJ*qsggE169lG+ax0YShN zh14<|W=l-m!M7uEfg))BFGvj$qp5URNOv*BJI=;X0u7&cRbWcxSnKr50o{Q=>uNE$ z<4RXYt)VUe5qewz@`_vq%8Fdu(M6r>h`Kx`OZR{);1O6TbD#vcdOt2?ZAR_jll;4E zi+lO^RQ_J5J=WtQ^nu+m#8_V5QN%s!WA|QM60G%h3j127s_3cp&4Opx`fFMS+uaGf0Yb&GC3r415vVJSnFaI7BrA z0ld|;=|&{ze63a_MC<)MA0;>?7-z}mBo$Ls+(uC^T8)O?Mn<%A9H!%GqLB~!e$QZV z0L3)=FO>c%4>zUk*^vagC=HZhT00Bp98eyp08|9oA)sNPOH6N< zfJT7I%xI4SjRB2=+yu}h&=k-#0}YJU^!5xhUieM@ve2Flh#)L|q_yW*radp_kLd|v z*&V%;tVq|A4s51_{zTdARpFCHM_8Y$%kN02#Q>9_LPg|CqSi`y{TD4J*u_ zB`~!F=1+T=mj=WvDa=2J%5k3MPC!2!Vdbab$z8NqPv{?2vY*Ce$6~TZLUtTvHC~Fz zPQ+wW3E4>{yBd?7ipi!EveQa-J0?35lg%V#XO-;JFcR%-;@%wG)A?MOZ6xi5$BOwNVI$MMc2sh~l$u zIpT+0nTJmR_73>XBSlZ%YbU|mlIVv6st;uS9(@nJHIMB)cs9gFT1W#TM5+R1Zv;+ zi%}oMP+w7~&l0GY5~zO}#8&rt4D*JH>eU2hDTbLE%zyW{G0X*pxu3usNnrAU?EGE~ z^HqiU=LBXsf%&IFFr(_}=)1p2RPL7mKg~PQ{qw6U>>{=ZeU|2nF#nc5Tf42BqHtUY zZST)d>?~Hn&QiQ(T-6NhDDdKHwtnCSI)5Tz+gFwB z55eZ~1pjM_9{_(W!T-AAcjwOWCH4l>PP7-!y(x%1%mMa|6E(rYsP&)dWOhGFAm0My zEO(>mo}YKf>hue;iu>bvS&br@LTxO4)0yO@O?8)T^J%I(jp_WA}8{qAc|s042w&mBt}G8jEXTaE+)jJm=e=sW@`rC z%hn#$x^H7dl5qLqdvD!eC`e<|3btjs>DznElWA`+Xzm3v)3k!X=B~`T(TZAx9@T8 zar4P`-LcK4#r?WlGlOm0RO03qvTvE!e6#MGt{0e=xojZ7S%`wO!*8`GoK?us_kH zw)~o;b#U9Y;=044s)9go>7cleT0_s_@o<96ywA@Z=X?-zT?!~ppY*gU?_n5sM*TYkA#nW(*I&{k^7T4*& z;lneqI-?sCuIdhdgQ>bJwkiCMfp^ic_EeG4lJ3fS&fFp3vlC94T5(E@&_l1?RtB)?NGOG*6RDV+cN>j984n9sk=64 zKgEu3tqlZy)*l51(aHJS-jfYAb8`W=Tpt{@2aFQNC3>drx#qrIg`pb_$6f-4z%dvY z8{t?xQ$QCEu2_ek*YLlFF!V~N>(3UqYBFW-|JZhU5$->yuRN!VdA`P!;!&X+ce0$4tn9Sm= z#I(wkAMkn;PEDsC$kdkOS%J*IhkJ@jv`jTQ^pVv4Lth&90SM_$uPM`36CjyPK?-ZO z)E#hq4!o5aZZ{pP3dmK*^X-KzGN+8Z4||!^_-(24hY$G(iqiKbeJtrZpTetrxbvc% z8V0KSRZUBS!N^zbg#Vs|9l+;t1%J{Iei6A7BN4(bHSOOi0jS@(4 zJ!QBKbtAXGj?0|49#%R+)X)!&-c;zQI}L>TulMwa)*~MYQA>0Jjv+$$Y!NVrfscM* uR2JbJ@-RxMrKs`y5c36%FT+W#V@bmI5Uw(nU2S-5&#~7qLq2XZS^YmfTDq_R diff --git a/src/__pycache__/read_input.cpython-39.pyc b/src/__pycache__/read_input.cpython-39.pyc deleted file mode 100644 index ac029ebfe32c42453c462dd5203890d252c7a953..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 5068 zcmcIoO>7)V74GWp`Ty~F?2Lb7cVL~3!H!K>76@MNDsjTD5UWUu0b8>^)vOR+ao}V)tPnzCytlosU{8BMdq5mHq7vc&65uQ36d}QT)$W<@#16BFozbhV ze)al$uijtvj9DmT75MD_+SsOZit-n{iGEbLc@4-Lhk+`T*vh6tRjR$KY!XWDDvSCF ztSE0*30_V=m*6@LbQ8$CpMfVy7o@vtUFm8-I#2>A36!F0Hx1VeP!=c$ln2=YP!VW^ zYTXjhC{UT|-7%nXpbH>Z0h$1s1e&TVH1R&^PE-BNZ}Jzp?u@T+g{BVG?kr7rFY=2$ za;(tIu6C4Wf9JD3^;qd?$BJCh;X23X>dLOlKG}PZ&(Z9Zg7SRkKw*ETIX*Mstto>! z&F`v5guT=$`3c|_rj$e2{dW{+Nmu;j5oxay1;2iZ7K8N)NWgu%mpKOh1I1BSm8eEK zA`yQSLsTP(KkOrxhD02dh#&H~UY2H$K|3wev3Fq?8)&T-(>^X`{}hqE5RuhmvK5e3 zSvDd&5s^*AWGAKU%MsbBh-@+@J1u3mBeF9Q*;Gt+R?5B;c%qw*tzCpQja7nJ%Q!>4 zWfsiIEch$W^dJM#Baj_32mHKW;B$xMkZdS#t#{}7Je{v7G%L0IgJ-A#d6DxWmV6## zjn6vekRNc7cckw9F4PP2pT+FDB<*Sq_SY%X&qa2ag=k-q`J=V%ezZT&FaFouc>$0U z?2}M4?CR+qpY9X<70Gb$(r=%W;|qYQv*}ny7Wx^fx5xW7UOqn~kBIa0fmZsW#B^eq z#R%rkkeIh4m{%m`Z(^7u5zLQ=!E_!+5ErFC%Q3`KA2HD`4u$xB1n~A{V`_um!#c4rMd1{ zjDHRIDF(eS?1i1nQZH1)A&FRp-CSB#*SgbSNl=f;_P_d;d^xn_IlcTPscmVXjeaFk z=aracFH6}Cu=qlZ|5eHV9`MIw{I5y=M}Y@=BlPQ3Jyy@`a+Bb7sC#<-SoW)@pEyeG zy%R(J1|X-|z0kMk#~dO$bw)({zB?nLfgjWDjlLh(&h}%d9jyC(%(_?LjD9a%joSwN z)Be_3;uFxP(URvx_B6{voBq(Zsh?%>*oIn<(t9@8G}|j3!x_$J;H z?>&YJSm#E(W5aWgIpHeTxXu$i$x}SdGd#<4JkJZf$VYgIkMc4f&08B7TiJui?}Uf5$SUxKTs)bz{LZ8lK^}zF`^%rrn?hnrtwW zHugQ_>{~2>aCiSlxM1%s*S8$MVcIcdk*osOVzJj)bfN}y7sxBYP>0vYx6)n0$xqZh z{jTzX!p4C8fhyFitAaG-ogALLva;hgtrgpH>!t^LrYjq6b-xKpyp@|SwN`Gjw&$Dn z^0#l@T=7_Sg;^#I{CBn07TGns;rUg!*;-!%4gUmEI&%@3HQTM4wzqzDAh?-)v}HbO zc%N<$1DhGX<40~M@OBc0S8*~0!)CJm*=V>i(9rp8G~5a`h;*3H|93yV9_nC8oTXq8 z8dSy!Ceq2j;lt7}dqyu-ip((#fXuOcm%vXIc;yW8Oa^6<2(97x=SJxPJS?Ns`PNVc zf(7X8H@7WjxHTwC*$DbI7*ybbx0$kGrBSsQ{ESnEgx)r3P_RY~>K+xHf?)NV&`m(+ zg2(WuC(FSAy~b|sn9!sw-xP_u*=(9D4+|&gn$EqU_Jf_na@_s;PN3T}>Wu@->6?J< z3_6ZsHyjJJpI~Ko)`xE+=pATw`P21Iv#A)8q!PXm++Z}>lg=C6`E_< z1p5Ifz{a34Yz6J>ndEtHgNju1dBy!pa6?es!Cq?2tc!$o@V4c!0<1kDFFqkltcczk zu8}|Mt^RRggnrANZVI(^uH1lp3C@D+Mzi7d&xc*B{m^C9SYF1{VR?CRROoKYg5#=c zdREP~sYqi_#E^=F=d(r&E-ky^i$u+KOQ|?Frp_1VUQNZHc7W z0!S(n5Wc!4G#fS!e}yO_#jKWXRsp$cyPmapNo1svx8b}JD!U;xcIOT&qbT{VAP)rD zWcaU#$Ofk?%E=}mt6t`_&?&TaSxSTKXY=q7GjPB^PM|?SPThdBk1h0XZ(~%@X>T}p zf+IX5ZzMkt9qpUp;VdJ?UVtYzv(+7|y88y6$%52`j&9miZI4wDm-vzFz+Nx1gYyOJ z&THZ6t?J~{e$(o#npJhRNb)2NUzL>LTb9@28S|`MJ%20%Mww(`ERiCaR4Xvzy$U1h zKx~viis~uDJ?ciC(Roznym?XT2vCE4sN|aDq3$FQ>c8CAA6SokBt$Kd3XB0l@N6EW y(H7hXdSyvIjXd-cYRRh%50OQct-?jlBSFG<4JIR%UTeB^-?r8B|gaUj{OZI3sg1k2sJ8Ifpna=j8(8oGi*k#Cf?Smk}3a zNuEJGCs*Vu;(56y&mvxs>+&4pqC793LA)r%SFFl~Cs>iSRi_n>{yZe;+rSrXA<$Uy zSxXBF*Q}w9yv=!p$UDe8oJWYf8(T+0V>io|b|_r6v}@vH%=wYKr36NO!U0WOu@k$~ zOWdfh;)<7uS|T1KVlNT7SYA=!XRCRn>pqpFe^i^r z$9qqZ(Aam$(ss$BT}gq^*3Pkz_A>U$Id>-KjdIc+wdKCHM!PhHHYG@f!l6$ya7!-L zrWdzuKzdx~=;=DlKALyZZ_VPdx#c;FdZI6nJt;3)RMT*V8J)q|;~X69P`MP)VN_g2 zB2pxf#Hn|q#EF`nN}w`SN}@?bJ@IO--cBuXAB1r&@fx*uyH@pZGz}cvm-wtqFh}hd z8ex3D9=2Nd+qGCV_lhWt=sYd^fqP!OD%OAOdTIsvdY($tzl(?|&&O~J3v5`!e2inz zj&`+&jb#q~VL@AOTKB#;oYQlpn!|a-P|4u}VyNPo!wTNC9*cYaa3yy1iWIAssXMf>gMS^a#@>;y^Lphh z!H4s^So5(MuBD}`Go@$snm&tFucJ1j*J)l-R-~|NO2Hi0Nlqq>an5OK8|d|8tYo;Z z&;4I}qV;*S&FTXB+W(IkpV7}qZ`V>k(yRK+ZU*rgNNSa(^P~Q!aZaDXnHUKnMbf#q zHe7)`^ZE=)L9eVq0>XN1_if1Xg3h!2J{dl%FJPo+Y1as`V}#f%LTH4Mf9&AI0+#Hj z0JBMZ8xtN?6c~*XH+)oUC4Q~f3p+CLRISqp6Q>n+7$O0JO4araQhKU8KpP(j*B+ox z)KzgP6Q|ehRW2skM%W2eE$*snA=+=hJ<3f$RV`u-YKekn3Q7p7cIORzzk+Wg@$22L zl2J9ct!nkJb~kD&P8K}PQrNbpxIy;KaY8$76MDS&XOGn{vj^V z&p9vfv)3J0tcVpKAL=iOd2t!Fm+S!T=Iy{KiWS>+a^fOteik3Kcp2ZdG*9tm$3^QE z^n1Zo*D&7Nzm;X6H0Q;IBbJc#Gx{V=(fWeC`vU~2NMj{7eV#Z z=g_B`F$IyT&@gObFH!0W4d|Q7lWLydZgwQ}oJKl8yNK3d*)GT;>G~S$`4Vz#;)1=>+h%(&{IOX9;hzHS;;_z@FAmNAu94ctI1O+9Us0nw$ph}aHKzqb= zP#R2?W~g)(l<{e;F-GRpsF|7mz8RJKpBfjlJ*7294Ni{Alr^ifsck}0n58}gIVXM# zcG+qqZl|l-iC0(MC<1D$Imvx>~>OC`%ScoE+DXMV4hnLMPOr26Y-~sdi7bzGg@fedLqckbv2@dpt^$O>26EsY<6hJ&sj98+o)n4CA_@Nt)e#x41uMA-^rnl}lC9oG$;{5_ z(lJKxo5-93!ArC^*aKEP0>VSc10V$`BwM1I5j4WL+vZw1vk`_6GGdhS9(ehF)Y}PF zXf|}>Xmj;JeBadnCHkdcP6)6@7sg*jS)JBC`XgsOzNzfNcVHj500Yd^Q+#aSy9rPh z8URa{$#DjEnGBsVYRki5w7+3CW>670d46P z;W569jDaWEw~rf5_7HR8J`>>XW$Wx{f3!;M+}DeA=fbMi zaJMQRL2K!!=L|o8Ec!3<&F*i91@y*RMvjj5_gJWhV05CN8_9$2cB}CDyEf`JV^0@j z{|MLwmzRgjE7BM&&AHS2f}GNwV_(b26=eatg=Xzy4!FT(UEIx3YAneMLnF!MROXN* z>c3C2#Hzlp5t@FoaD&#l$!jfQwN*@e+Awdm z``7!I;rd;)SsE7PJhYzHH;0k!{=1O2dkFusZA1STw$v{`87GdHee4r*NndUJ+s+Sv z|1bXY)7Ng9+x{<6!ktqRB|fWjX&SHrD~0;-{d$Eog|aiSA$W@iEqX;SiPEUcJ7$b)Yw&5{aU>q_CR>U zy?WKHVO(_y^Cw~#3|BR4t&>PJ$4(qnxj~OZe}lFW5fZp@dx-&g9>82gE?{s9go{%8 z2(EpxB+39RBCti9h$@j#AKwz{%b>TWrDf3E<@0o|yj77kjx(!U5LCdfl6)D9dEKxnCpa@-PNeh~!L zf+BcZY7;Z5FHo?CplWZlKxSvpw&I%ubkv`oL8WY@@^Src)klj1cTsmDz% zFr0`G&|+p}vp{XIcZk+jUq<;1<`dM8*-;67VfFTh#^)>bDS~_8~vGbriP|@WZt!};6ie5uv)WeW<;_~L1aGI%4s90tj zH=sD|B;NbvT@!D+)y13N#03{v;c#`1yX`2M8_&s{o%M|k+io?hb_-Y5CgO*P>rH$Q znYd5(s?77^iPaZrZ)Gghg5}QJbpOJ&6n4J^TQ1rK#4kWSmxTH)w3$)P916yk`dwt8 zWS~m0HaImrDUm`tq>ylRc({X`8QJUaOxx?wLu*ndkCh8<6bcAp&Xy~ndKuJ#j;D$h z!~y41#TjOa6vK&sJI+cEPo2nZJ1dw;SaI$se+-M$&Jo-h(RUr3GX&@b?XcqTeakoi zI$%Dx0fFCJN_AEF)pFeJmfM5+PI)Z+sJx9DI+Akd3I?_(zC{E!BkB3um z^d3lHN?yk>K_CS42ML7RCK#eSHlX4W?&NfzPbn?Os7LG_1`HE?S<~J*JX8$r{-u72 zfdcRM?8gZm5|KMie7}YY_;i+Mof
=`kVrQtEC`Ot{(9=y`S@@}1$h4eS5%N-klGtM za|Ai!?E_a=iT_m`!|p(6*HG`G{=IQM|G#Xe9u}wcR6SbE)C2C-cgFQR!%Te!^?yB8 zf3kmo`hOhP)BLmJXHma4RX?+Se(X;X1H||iT#N@&N^T&l{aqiV@p%wMT-A?ILVIS^ zK^mF>`kg-c=V1MWODJU*;vy%;IS;aPwF3eU?c_1&UkoFh`UgC#=FMEu2WAYJV>gz^ z97ydtSjC@W6(6K?Q~6}>g3jQqJ{*#>%;ASC{h=-MOvyh%Sz)S-DDfvIkBbTIDbMSr z)z1;dMo3qTIikWdYoYzl<7&$fx1hpP`(xVw0oL~kt&ye0{|?$F^jp}fIvdyajISbA zn;^SSWd~Jf?>166cUyk{%nTU_&-*>fCT^m diff --git a/src/__pycache__/utils.cpython-39.pyc b/src/__pycache__/utils.cpython-39.pyc deleted file mode 100644 index 9300972f0696a318a254f573a9a6db339cb353ae..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 8202 zcmb7JO>7*=b?&PEnVvs#I2;ZsN~@h!>_8NDrCoWm-YAN!HDwuAcDaJ=7{(y7Fq&!( zN1Wg49VAx2$ z_qt~|Ln(2>LD%1_-&e1`_f@e`EaoiyeQwiz_~SP$>t~di{gcSt#25ZK0%a-DvFd_j zS(hB!bz4bge`wVm?vPm5oc6R*_4xk~*heLcFBJSFMfn&oCovt4TAQ{JBWbw}~&@LSQiC zGnNq)u32M=yyQGW#!2F1jQP~r(gM9c?_m0SJhs2lX6@XBUi z1y?rpei$`6uf6xB%`2g9UKvI0PI!5+AN#5eL9ZFmqCy{m<(33~QZFFg^r$5LBWV^N zuRTV>VBHnV$cja)QUakg_OVcM8EfU7GnaEGIVC4`xvey)m!?ps0+oDUSCMb#sWkn&fSB^U4Yx7DmNm{t*al-8 z$GBKn=D-^ljCISp|IxT;ie#GO1;jAP@gicF;<#i4JfvhqkLq3&NOdiI9)~13PVf+@h zRanDo!*#4@C9)&OrmYH#?{C+RV7xfM0Oytf}f99BlpNNc~d<@@L~T6 z=6o#1Yia4~TyT#5I*m)liWF8&DICXj(o>Q~KW7cq_09S*W-?wkXaBD) zQTrU~W=#QY<^N;Em&{Ab-L>>D&8k`16Np!!sa4j_FNdE-IkSR2Ni>8MY3KghxC(vd z%?fG3RM(&ZVI9h$gf7pUJnQd=)bw>`_A0llePZ3dC!k|y3O-rJ~_09w`b=3>fD@X_4|C=f53b{LucKb$Bw^?-F$8* z0bjM$=4!CHSNM#*!u$OS@AoSw`>nhdfD=|?5C6mD0eP7j$a}`%N4VGe!^-!*{7+9l z{6Z`?HwuhKu@meyIm5u*;ddCbR+8PxDf7lzdXsEfT~`?9CVq2WeO??;1J)z_eFf0vDfVPwF>LG zZQW>ot-smOQEa1IO`b{<_0R;CbsPQ1Vib$LSbR`-qyDX+6N`2%9>r3>f$wVF)^A*m zrE1`Vi&gb(FZf`Sdyi!!iiNH_`tkmg(=hHoY_}N3&taUfgTNAwU6#J!AL0W2>~jJ? zxo$h6Dykkn$S;WnaT(G}(nq}o>Dwhym5!Yg7a;jfd{pBtd{@&v#g}adwX0}%(b2z! z@#ZKBZ5??j_78WEV5Gl^(=3&6Eu43zET6Y*>pHxKl4%)E=d{dD%W&S3vY}%@Pc$(S zUjK%zFCiW0S^e!^BMQP~ugHPTIJs;D^9VbT4ug}!Y(8d0E=^XYVh2AE4d*91Wf z^t)hy`Yp7nXOemFjHcWhTs^)s&Z1| z8>5+w@FnN6m zx{kdLj%%kCJH5W{#%@#h!w_Jun@nm|Khle+riqPi^dOP{J=&!RPv^f%x5)QU)})zF{whdJuP+@q6s!WbVjz6F zhlle1n?P$}fW34sInLmMlc61kv@-US)s3?e`2k?P0figqm>kiy=R3Im6d3J2B0sxQ z-<#Gfpe9|4T;Ta1;llBAQ#3iy?R(0@CCNtX_3^Z^?JD8lgdBvMC)oc-_Fy6ceg`D25mKlJ~=L0c{J3CJx|>cAZ0^zX@t5hA2<;K$pB# zIR9M<+0Do`rN}!18o~b+;r~k118cGP{EDEVbj?`NGIAxgfOW!}b}@#&iA!8*Pf%*A z!}F;Qms4#+hme0lTEwhAF$j}(YT-J~^Cr)=gtb<%cX^O~qijpq{L+?w8wSXRugTx? zHUAkop|7?6ZRcly`gi~N>5bcoza^{Hbj}i={X3NUT?Efa5R3P5`9C@TbZ}x9zA$*4 zxZ)ZL^i2xhL9kKb;}rXY207wh2zEJ1#@Sv!`bzJ%ZnWxNo8wsY;!G-!{YlLezl@i$ zW z)G@20E3pkpr$ne?S1qv#|VvkWfeoAt@Eg_$JEd@ScEn%8V*#tMIb+VdC%7 z#p5ABvNr=CI&zMg5hYw%*~!Dn{g`+Uy*1+5Bw^9SmOo~t;RF+IK|zTEDea*@M&Dg< zCgQ*2Oj_biOZ>+(_`W|D zE$c7k&;cQu`Y(D~z5M!%o!a#dSaGx89o$$#@BfZ|tS7(0C#4s(TRRW?BOPA8-sv|R zo$v+{liJCvAp;A|YL1~ZMPLe+*~Wz^3VN~o0r}O~-R|`91bE_n3yHD0xEgoExHuh? zxjNe#A8OfY*JTGc*f!!vh?{MEkC?EZtWzhiji+W`!u-M-Eexl2LDF>$w^W?@GR|>H z77$;AX_kfFMV&drB($8K((fY!0|Q0ENt2Pl3lkZmO~wdkhj%ZsaIV=t&RYF%8|x}$Ebzm57EMONwwDc6K` zGlGkAmKVwB;qT3c62PA_uE*U<+gkOiBS-PP-h7-1C&T3jzE&>vKSrY#j^|HPT7L+fe&a zqP|DSp4Q<~Z`AAr4gCO(Mx7A1oOei%3Di4$w%C_KgXfXcd+1}%@6pMh>hDfMiftms z`0y+-_VD3q3gG?(S5qIa=u%%ZnIq_tpF%)^kMRR6uF+)>(r?U2ojEB^o8Fs|dUMhY zq<=jl_2;BnNPjRR&CN;kkp6T=T9}g-Azho1E}TgDwLirOkl`=Fp|~KUSPxFGw@1{m zGW><{j*zESGRh#0od?ZJU;J}8?Sw-pWq08qhl@d)bD=B8*r42{$iwK0e9fTZ@)v8nuI?1IT)Z$2E8W6a@)t2>DJ5)<){P*#{Jqt$$r zm%_Yutzy5j!@Z4=h5 zY}M_}tB=8G$oiSVazLf~C?!gJBBIhKlhUcEZa;pX>T(a$Zm5*%P763yxA(bp+BlID zB+d29@+k#xy|#JoU;p{7^)JkDZ=z)znS@SHo#zZn(xGLixv?;Fa(%q?c7sUo*Q;r& z{UlKQkej?y*2}Zy31Ze4N4@rN6nupUJbBrACei_3Emhq`?e_=z4h=?&sN0R`4JOr; zx*ChYX*ehmaTrCxgUL+p5gOb>Tf7X@54h^s^g2v$z%Trmd*_72lL(kK`C%aOXV8)u z_C>VHiPupI+PIHiH_EZw8+8Z!iATAKBHg2a8AX*+^i0gpt^5MV&no=Fz!!8T;!Hpj zzeM=XAjQ4UERTH$(+U527Y{+5;0BhK@QC!aJMS0#-}kfrihs_}BYwv(`Ro3wzx3bu CnwI7O diff --git a/src/read_input.py b/src/read_input.py index 9dd2f91..0e0bb37 100644 --- a/src/read_input.py +++ b/src/read_input.py @@ -9,9 +9,6 @@ def read_input(inputFile): X = [] Y = [] Z = [] - xOy = [] - xOz = [] - yOz = [] symmetry = [] outputFile = "" pattern = [] @@ -70,27 +67,6 @@ def read_input(inputFile): elif ls[0].casefold() == 'z_axis': ls.pop(0) Z = [i for i in ls] - elif ls[0].casefold() == "xoy": - line = f.readline() - ls = line.split() - xOy.append(ls) - line = f.readline() - ls = line.split() - xOy.append(ls) - elif ls[0].casefold() == "xoz": - line = f.readline() - ls = line.split() - xOz.append(ls) - line = f.readline() - ls = line.split() - xOz.append(ls) - elif ls[0].casefold() == "yoz": - line = f.readline() - ls = line.split() - yOz.append(ls) - line = f.readline() - ls = line.split() - yOz.append(ls) elif ls[0].casefold() == 'symmetry': ls.pop(0) symmetry = [i for i in ls] @@ -216,4 +192,4 @@ def read_input(inputFile): print("Bad input : missing the keyword -- %s --"%t) sys.exit() - return rB , rPP, center, X, Y, Z, xOy, yOz, xOz, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translate + return rB , rPP, center, X, Y, Z, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translate diff --git a/src/read_input.py~ b/src/read_input.py~ deleted file mode 100644 index 1682874..0000000 --- a/src/read_input.py~ +++ /dev/null @@ -1,196 +0,0 @@ -import sys - -# Parses the input file -def read_input(inputFile): - - rB = 0.0 - rPP = 0.0 - center = [] - X = [] - Y = [] - Z = [] - symmetry = [] - outputFile = "" - pattern = [] - npattern = [] - atoms = [] - dist = [] - a = 0.0 - b = 0.0 - c = 0.0 - alpha = 90.0 - beta = 90.0 - gamma = 90.0 - showBath = False - evjen = False - showFrag = False - notInPseudo = [] - notInFrag = [] - symGenerator = [] - generator = [] - translate = [0.0,0.0,0.0] - - checkInput = {'bath':False,'pseudo':False,'output':False,'pattern':False,'npattern':False,'a':False,'b':False,'c':False,'atoms':False,'symmetry_generator':False,'generator':False} - - f = open(inputFile,'r') - - line = 'x' - - while line.casefold() != 'end_of_input': - line = f.readline().strip() - ls = line.split() - print(ls) - if ls == []: - continue - if ls[0].casefold() in checkInput: - checkInput[ls[0].casefold()] = True - if ls[0].casefold() == 'bath': - try: - rB = float(ls[1]) - except ValueError: - print("Error while parsing the input file : %s is not a valid bath radius value"%(ls[1])) - sys.exit() - elif ls[0].casefold() == 'pseudo': - try: - rPP = float(ls[1]) - except ValueError: - print("Error while parsing the input file : %s is not a valid pseudopotential radius value"%(ls[1])) - sys.exit() - elif ls[0].casefold() == 'center': - ls.pop(0) - center = [i for i in ls] - elif ls[0].casefold() == 'x_axis': - ls.pop(0) - X = [i for i in ls] - elif ls[0].casefold() == 'y_axis': - ls.pop(0) - Y = [i for i in ls] - elif ls[0].casefold() == 'z_axis': - ls.pop(0) - Z = [i for i in ls] - elif ls[0].casefold() == 'symmetry': - ls.pop(0) - symmetry = [i for i in ls] - elif ls[0].casefold() == 'output': - outputFile = ls[1] - elif ls[0].casefold() == 'pattern': - line = f.readline() - while line.strip().casefold() != 'end': - pattern.append([]) - ls = line.split() - for i in range(len(ls)): - if i%2 == 0: - pattern[-1].append(int(ls[i])) - else: - pattern[-1].append(ls[i]) - line = f.readline() - elif ls[0].casefold() == 'npattern': - ls.pop(0) - try: - npattern = [int(i) for i in ls] - except ValueError: - print("Error while parsing the input file : the number of patterns is not valid %s"%(line)) - sys.exit() - elif ls[0].casefold() == 'lattice': - line = f.readline() - while line.strip().casefold() != 'end': - ls = line.split() - if ls[0].casefold() in checkInput: - checkInput[ls[0].casefold()] = True - if ls[0].casefold() == 'a': - try: - a = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - elif ls[0].casefold() == 'b': - try: - b = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - elif ls[0].casefold() == 'c': - try: - c = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - elif ls[0].casefold() == 'alpha': - try: - alpha = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - elif ls[0].casefold() == 'beta': - try: - beta = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - elif ls[0].casefold() == 'gamma': - try: - gamma = float(ls[1]) - except ValueError: - print("Error while parsing the input file : bad value for the lattice parameter %s"%ls[1]) - sys.exit() - line = f.readline() - elif ls[0].casefold() == 'atoms': - line = f.readline() - while line.strip().casefold() != 'end': - ls = line.split() - if(len(ls)) != 4: - print("Error while parsing the input file : not enough values given for the atom in line %s"%line) - sys.exit() - try: - atoms.append([ls[0], float(ls[1]), int(ls[2]), float(ls[3])]) - except ValueError: - print("Error while parsing the input file : bad value for the atom %s"%line) - line = f.readline() - elif ls[0].casefold() == 'show_bath': - showBath = True - elif ls[0].casefold() == 'translate': - ls.pop(0) - try: - translate = [float(ls[0]),float(ls[1]),float(ls[2])] - except ValueError: - print("Error while parsing the input file : the translation vector is not valid %s"%line) - sys.exit() - elif ls[0].casefold() == 'not_in_pseudo': - ls.pop(0) - notInPseudo = [i for i in ls] - elif ls[0].casefold() == 'show_frag': - showFrag = True - elif ls[0].casefold() == 'evjen': - evjen = True - elif ls[0].casefold() == 'symmetry_generator': - line = f.readline().replace("'","") - while line.strip().casefold() != 'end': - symGenerator.append(line.split(',')) - line = f.readline().replace("'","") - elif ls[0].casefold() == 'generator': - line = f.readline() - while line.strip().casefold() != 'end': - ls = line.split() - try: - generator.append([ls[0],float(ls[1]),float(ls[2]),float(ls[3])]) - except ValueError: - print("Error while parsing the input file : bad value for the generator atom %s"%line) - sys.exit() - line = f.readline() - elif ls[0].casefold() == 'not_in_frag': - line = f.readline() - while line.strip().casefold() != 'end': - ls = line.split() - try: - notInFrag.append([float(ls[0]),float(ls[1]),float(ls[2])]) - except ValueError: - print("Error while parsing the input file : bad value for the atom %s"%line) - sys.exit() - line = f.readline() - f.close() - for t in checkInput: - if checkInput[t] == False: - print("Bad input : missing the keyword -- %s --"%t) - sys.exit() - - return rB , rPP, center, X, Y, Z, symmetry, outputFile, pattern, npattern , atoms, dist, a, b, c, alpha, beta, gamma, showBath, evjen, showFrag, notInPseudo, notInFrag, symGenerator, generator, translate diff --git a/src/utils.py b/src/utils.py index 0a89a33..512f712 100644 --- a/src/utils.py +++ b/src/utils.py @@ -119,7 +119,7 @@ def translate(v,coordinates): # Finds the point at the center of the given atoms that are the # closest to the origin -def find_center(centerList, coordinates, without=[]): +def find_center(centerList, coordinates): centers = [] for i in range(len(centerList)): @@ -129,14 +129,6 @@ def find_center(centerList, coordinates, without=[]): c.append(distance(c,[0,0,0])) # Computing the distance to the origin for at in coordinates: - w = True - for i in without: - d = distance(at, i) - if d < 1e-6: - w = False - break - if not w: - continue if at[3] in centerList: centers = sorted(centers, key=operator.itemgetter(3)) # Sorting the list with respect to the distance to the origin d = distance(at,[0,0,0]) diff --git a/src/utils.py~ b/src/utils.py~ deleted file mode 100644 index 348e604..0000000 --- a/src/utils.py~ +++ /dev/null @@ -1,390 +0,0 @@ -import numpy as np -import operator -import sys - -def distance(a,b): - # Returns the 3D distance between a and b where - # a and b are array where x, y and z are at the - # position 0, 1 and 2 - - x = a[0]-b[0] - y = a[1]-b[1] - z = a[2]-b[2] - - return np.sqrt(x**2+y**2+z**2) - -def get_cell_matrix(a,b,c,alpha,beta,gamma): - # Computing the volume of the primitive cell - omega = a*b*c*np.sqrt(1-np.cos(alpha)**2-np.cos(beta)**2-np.cos(gamma)**2+2*np.cos(alpha)*np.cos(beta)*np.cos(gamma)) - - # Computing the matrix - M = [ - [a,b*np.cos(gamma),c*np.cos(beta)], - [0,b*np.sin(gamma),c*(np.cos(alpha)-np.cos(beta)*np.cos(gamma))/(np.sin(gamma))], - [0,0,omega/(a*b*np.sin(gamma))] - ] - return M - -def big_cell(generator,symGenerator,a,b,c,alpha,beta,gamma,nA,nB,nC): - coords = [] - - # Computing the matrix converting fractional to cartesian - fracToCart = get_cell_matrix(a,b,c,alpha,beta,gamma) - - for gen in generator: - x = gen[1] - y = gen[2] - z = gen[3] - - for sym in symGenerator: - u = eval(sym[0]) - v = eval(sym[1]) - w = eval(sym[2]) - - # Making sure the value is within the range [0,1] - u = u + 1*(u<0) - 1*(u>1) - v = v + 1*(v<0) - 1*(v>1) - w = w + 1*(w<0) - 1*(w>1) - - coords.append([u,v,w,gen[0]]) - - # Deleting the redundant atoms - toDel = [] - for i in range(len(coords)-1): - for j in range(i+1,len(coords)): - # Computing the distance using the minimum image convention - # as described in Appendix B equation 9 of - # "Statistical Mechanics : Theory and Molecular Simulations - # Mark E. Tuckerman" - - r1 = np.array(coords[i][:3]) - r2 = np.array(coords[j][:3]) - r12 = r1-r2 - da = np.sqrt(r12[0]**2+r12[1]**2+r12[2]**2) - r12 = r12 - np.round(r12) - db = da - np.sqrt(r12[0]**2+r12[1]**2+r12[2]**2) - r12 = np.matmul(fracToCart,r12) - d = np.sqrt(r12[0]**2+r12[1]**2+r12[2]**2) - - if(d<1e-2): - # We check if we don't already want to delete this atom - if j not in toDel: - toDel.append(j) - - toDel = sorted(toDel) - - # We delete the atoms in the list - for i in range(len(toDel)): - coords.pop(toDel[i]-i) - - newCoords = [] - - # We replicate the cell nA, nB, nC times - for at in coords: - newCoords.append([at[0],at[1],at[2],at[3]]) - for a in range(1,nA): - newCoords.append([at[0]+a,at[1],at[2],at[3]]) - for b in range(1,nB): - newCoords.append([at[0]+a,at[1]+b,at[2],at[3]]) - for c in range(1,nC): - newCoords.append([at[0]+a,at[1]+b,at[2]+c,at[3]]) - for c in range(1,nC): - newCoords.append([at[0]+a,at[1],at[2]+c,at[3]]) - for b in range(1,nB): - newCoords.append([at[0],at[1]+b,at[2],at[3]]) - for c in range(1,nC): - newCoords.append([at[0],at[1]+b,at[2]+c,at[3]]) - for c in range(1,nC): - newCoords.append([at[0],at[1],at[2]+c,at[3]]) - - # Now we convert the fractionnal coordinates to cartesian coordinates - coords = [] - - for at in newCoords: - r = [at[0],at[1],at[2]] - rxyz = np.matmul(fracToCart,r) - coords.append([rxyz[0],rxyz[1],rxyz[2],at[3],'C']) - - - # Returns the list of the atoms [x,y,z,label,second_label] - return coords - -# Translates all the coordinates with the vector v -def translate(v,coordinates): - for c in coordinates: - c[0] += v[0] - c[1] += v[1] - c[2] += v[2] - return coordinates - -# Finds the point at the center of the given atoms that are the -# closest to the origin -def find_center(centerList, coordinates): - - centers = [] - for i in range(len(centerList)): - centers.append([100,100,100]) # Setting a large value for each center - - for c in centers: - c.append(distance(c,[0,0,0])) # Computing the distance to the origin - - for at in coordinates: - if at[3] in centerList: - centers = sorted(centers, key=operator.itemgetter(3)) # Sorting the list with respect to the distance to the origin - d = distance(at,[0,0,0]) - if d <= centers[-1][-1] and d > 0.0: - centers[-1] = [at[0],at[1],at[2],d] - - center = np.mean(centers,axis=0)[:3] # Computing the barycenter - - return center - -# Defines a rotation matrix that will put r1 at the position r2 -def rotation_matrix(r1,r2): - - r1 = np.array(r1)/np.linalg.norm(r1) - r2 = np.array(r2)/np.linalg.norm(r2) - - # Computing the cross product which is the vector around which - # the rotation is done - crossProduct = np.cross(r1,r2) - crossProduct = crossProduct/np.linalg.norm(crossProduct) - - # Computing the angle of the rotation - a = np.arccos(np.dot(r1,r2)) - - c = np.cos(a) - s = np.sin(a) - x = crossProduct[0] - y = crossProduct[1] - z = crossProduct[2] - - M = [ - [x**2*(1-c)+c,x*y*(1-c)-z*s,x*z*(1-c)+y*s], - [x*y*(1-c)+z*s,y**2*(1-c)+c,y*z*(1-c)-x*s], - [x*z*(1-c)-y*s,y*z*(1-c)+x*s,z**2*(1-c)+c] - ] - - return M - -# Rotates all the coordinates using the rotation matric M -def rotate(M,coordinates): - for i in range(len(coordinates)): - r = [coordinates[i][0],coordinates[i][1],coordinates[i][2]] - rV = np.matmul(M,r) - coordinates[i][0] = rV[0] - coordinates[i][1] = rV[1] - coordinates[i][2] = rV[2] - - return coordinates - -# Cuts a sphere centered on the origin in the coordinates -def cut_sphere(coordinates,r): - sphere = [] - for i in range(len(coordinates)): - if distance(coordinates[i],[0,0,0]) <= r: - sphere.append(coordinates[i]) - - return sphere - -# Finds the fragment in the coordinates -def find_fragment(coordinates, patterns, npatterns,notInFrag): - - inFrag = [] - - for n in range(len(patterns)): - pattern = patterns[n] - npattern = npatterns[n] - for i in range(npattern): - c = [100,100,100] - dc = distance([0,0,0],c) - - inPattern = [] - # Finding the closest atom of the first type in the pattern - for at in coordinates: - if at[3] == pattern[1]: - d = distance([0,0,0],at) - if d > 10: - break - if d < dc : - accept = True - for exc in notInFrag: - d = distance(exc,at) - if d < 1e-5: - accept = False - if accept and coordinates.index(at) not in inFrag: - c = [at[0],at[1],at[2],0.0, coordinates.index(at)] - dc = distance([0,0,0],c) - # Finding the rest of the pattern around the atom previously found - atIn = [] - for j in range(0,len(pattern),2): - d = distance(c,[100,100,100]) - # Initializing the atoms - for k in range(pattern[j]): - atIn.append([100,100,100,d]) - - for at in coordinates: - if distance(at,[0,0,0]) > 10: - break - if at[3] == pattern[j+1]: - atIn = sorted(atIn,key=operator.itemgetter(3)) - d = distance(at,c) - trial = [at[0],at[1],at[2],d,coordinates.index(at)] - if d < atIn[-1][3] and trial not in atIn: - accept = True - for exc in notInFrag: - d = distance(exc,trial) - if d < 1e-5: - accept = False - if accept: - atIn[-1] = trial - for at in atIn: - inPattern.append(at[4]) - - for at in inPattern: - if at not in inFrag: - inFrag.append(at) - - for at in inFrag: - coordinates[at][4] = 'O' - - return len(inFrag), coordinates - - -# Finds the pseudopotential layer around -# the fragment -def find_pseudo(coordinates, rPP, notInPseudo): - - for at in coordinates: - if at[4] != 'O': - continue - for i in range(len(coordinates)): - if coordinates[i][4] != 'C': - continue - d = distance(at,coordinates[i]) - if d < rPP: - coordinates[i][4] = 'Cl' - - return coordinates - -# Creates lists containing the neighbours of each -# atom -def find_neighbours(coordinates, atoms): - neighbourList = [[] for i in range(len(coordinates))] - - atoms = np.array(atoms).flatten() - - for i in range(len(coordinates)-1): - for j in range(i+1,len(coordinates)): - li = coordinates[i][3] # Label of the atom i - lj = coordinates[j][3] # Label of the atom j - - ii = np.where(atoms==li)[0] - jj = np.where(atoms==lj)[0] - - ci = float(atoms[ii+1]) # Charge of the atom i - cj = float(atoms[jj+1]) # Charge of the atom j - - if ci*cj < 0: # Checking if the charges have opposite signs - d = distance(coordinates[i],coordinates[j]) - - if d < float(atoms[ii+3]) and d < float(atoms[jj+3]): - neighbourList[i].append(j) - neighbourList[j].append(i) - return neighbourList - -# For each atom, finds if it has the correct number of neighbours, -# if not, modify its charge -def evjen_charges(coordinates,atoms): - neighbourList = find_neighbours(coordinates,atoms) - - atoms = np.array(atoms).flatten() - - charges = [] - - for i in range(len(coordinates)): - li = coordinates[i][3] - ii = np.where(atoms==li)[0] - - nr = len(neighbourList[i]) - nt = int(atoms[ii+2]) - ci = float(atoms[ii+1]) - - if nr > nt: - print("Error : too much neighbours for atom n°%i, count %i neighbours where it should have a maximum of %i"%(i,nr,nt)) - sys.exit() - charges.append(ci*nr/nt) - - return charges - -# Computes the nuclear repulsion -def nuclear_repulsion(coordinates,charges): - - rep = 0.0 - - for i in range(len(coordinates)-1): - for j in range(i+1,len(coordinates)): - rij = distance(coordinates[i],coordinates[j]) - ci = charges[i] - cj = charges[j] - - if(rij < 1): - print(i,j,"\n",coordinates[i],"\n",coordinates[j],"\n",rij) - - rep += (ci*cj)/rij - return rep - -# Computes the symmetry in the whole system -def compute_symmetry(coordinates,charges,symmetry): - symmetrizedCoordinates = [] - symmetrizedCharges = [] - uniqueIndexList = [] # The list containing the indexes of the unique atoms - - treated = [] # Will store the index of the atoms already treated - - symOp = [] - - # Storing all the symmetry operations - for s in symmetry: - if s == 'C2x': - symOp.append(np.array([1,-1,-1])) - elif s == 'C2y': - symOp.append(np.array([-1,1,-1])) - elif s == 'C2z': - symOp.append(np.array([-1,-1,1])) - elif s == 'xOy': - symOp.append(np.array([1,1,-1])) - elif s == 'xOz': - symOp.append(np.array([1,-1,1])) - elif s == 'yOz': - symOp.append(np.array([-1,1,1])) - elif s == 'i': - symOp.append(np.array([-1,-1,-1])) - - for i in range(len(coordinates)): - print(i) - if i in treated: - continue - - treated.append(i) - at1 = np.array(coordinates[i][:3]) - symmetrizedCoordinates.append(coordinates[i]) - symmetrizedCharges.append(charges[i]) - uniqueIndexList.append(len(symmetrizedCoordinates)-1) - - for j in range(len(coordinates)): - if j in treated or coordinates[i][3] != coordinates[j][3]: - continue - - at2 = np.array(coordinates[j][:3]) - - for s in symOp: - if distance(at2, at1*s) < 5: - if distance(at1,at1*s) > 1e-4 and distance(at2,at1*s) < 1e-4: # Checking if op.at1 != at1 and that op.at2 = at1 - p = at1*s - treated.append(j) - symmetrizedCoordinates.append([p[0],p[1],p[2],coordinates[i][3],coordinates[i][4]]) - symmetrizedCharges.append(charges[i]) - break - - return symmetrizedCoordinates,symmetrizedCharges,uniqueIndexList