From 7899cb4decfb8987cceb7785c789554bfc7b1c17 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Fran=C3=A7ois=20Villemot?= <francois.villemot@espci.fr>
Date: Thu, 15 Jan 2026 12:01:39 +0100
Subject: [PATCH 1/3] feat: script to run validation set

---
 validation/lj-mixture/calculated-energies.csv |  24 ++
 validation/lj-mixture/correlation-plot.jpeg   | Bin 0 -> 21850 bytes
 validation/lj-mixture/reference-data.csv      |  24 ++
 validation/lj-mixture/run-validation.py       | 224 ++++++++++++++++++
 4 files changed, 272 insertions(+)
 create mode 100644 validation/lj-mixture/calculated-energies.csv
 create mode 100644 validation/lj-mixture/correlation-plot.jpeg
 create mode 100644 validation/lj-mixture/reference-data.csv
 create mode 100644 validation/lj-mixture/run-validation.py

diff --git a/validation/lj-mixture/calculated-energies.csv b/validation/lj-mixture/calculated-energies.csv
new file mode 100644
index 0000000..05385a6
--- /dev/null
+++ b/validation/lj-mixture/calculated-energies.csv
@@ -0,0 +1,24 @@
+,t,x,density,u,energy,energy_err,lr_correction,acceptance_rate
+0,1.2183,0.5,0.8,-6.991912,-6.824332775509206,0.007332089000306406,-151.9222906615277,0.523858
+1,1.2183,0.5,0.7,-6.28399,-6.15696393793348,0.0115122047621387,-132.93200432883674,0.605872
+2,1.2183,0.5,0.65,-5.86756,-5.747027114561171,0.008533981399021708,-123.43686116249125,0.643768
+3,1.5096,0.25,0.8,-5.87326,-5.728584068606758,0.021466385837235483,-127.1877273361547,0.580572
+4,1.5096,0.25,0.7,-5.2830446,-5.1486481236643895,0.015199454304553882,-111.28926141913534,0.65148
+5,1.5096,0.25,0.6,-4.57915,-4.457067486522425,0.012901052999909882,-95.39079550211606,0.713526
+6,1.5096,0.25,0.5,-3.84601,-3.7550796383339944,0.016097072216878654,-79.4923295850967,0.767107
+7,1.5096,0.25,0.4,-3.145625,-3.031679704426388,0.01957993121178323,-63.59386366807736,0.812887
+8,1.5096,0.25,0.3,-2.446133,-2.29186093973726,0.019961432466160035,-47.69539775105803,0.853689
+9,1.5096,0.25,0.2,-1.71833,-1.5204237232517432,0.018705478102023994,-31.796931834038674,0.897772
+10,1.5096,0.25,0.1,-0.88341,-0.7731743940092469,0.012321114326961577,-15.898465917019347,0.94422
+11,1.5096,0.75,0.8,-7.528813,-7.387557821039766,0.019429067147871005,-178.9234450957788,0.515488
+12,1.5096,0.75,0.7,-6.89008,-6.705393398859641,0.024860213510585486,-156.55801445880638,0.598998
+13,1.5096,0.75,0.6,-6.0058145,-5.844619979911088,0.01211174105592481,-134.19258382183412,0.674251
+14,1.5096,0.75,0.55,-5.52616,-5.422084348147788,0.005698415488936215,-123.00986850334797,0.708598
+15,1.6438,0.5,0.8,-6.54401,-6.36678624267275,0.008473648040505592,-151.9222906615277,0.557943
+16,1.6438,0.5,0.7,-5.95657,-5.844425928860567,0.013831459037555685,-132.93200432883674,0.634995
+17,1.6438,0.5,0.6,-5.18836,-5.085670690270235,0.010455755569991681,-113.94171799614583,0.702152
+18,1.6438,0.5,0.5,-4.35859,-4.21187686179021,0.010561762660606716,-94.95143166345483,0.758437
+19,1.6438,0.5,0.4,-3.549019,-3.4182878823413763,0.011275107483718705,-75.96114533076386,0.807701
+20,1.6438,0.5,0.3,-2.768327,-2.575420052801317,0.010850764609473712,-56.97085899807291,0.853085
+21,1.6438,0.5,0.2,-1.936931,-1.7475874342217939,0.01920265179952064,-37.98057266538193,0.89591
+22,1.6438,0.5,0.1,-1.004111,-0.8560035870798322,0.012956204273870288,-18.99028633269097,0.942764
diff --git a/validation/lj-mixture/correlation-plot.jpeg b/validation/lj-mixture/correlation-plot.jpeg
new file mode 100644
index 0000000000000000000000000000000000000000..abfe4c1b4a743e8271f65a8a21c377189aa7d415
GIT binary patch
literal 21850
zcmeHv1yo#Fy6(XV5&{8&dxE<=30e@`AwYslfZ$Fc0fJkA-~@MqdvJG$!X-Ea2(GWv
zJ$I%%(|zAech6nxt#_-|URA4V*EzNK`Stt1|J+UAEddXtB|s7Y3=9ll2z>%~vw$do
zjEIPYh=7cQgmmv7G71_tIvOe}8a~!T3~XWoQW9bUA|f&>CTcPY21+6#8g5z!7FG@p
z4pM4fK^`^%CUy?C?+1aockdn=DjFU-IvyK25joqx{kUreFp*)FV7=gAC;(VY7&uIr
zyLNyS0ALWH|MvSI{@Viv77iW(5efMo3MzC$<pTf~1`ZAu9u5Hk9v(W|1NuDxkBNXq
z&iWMbq2enf3L9)T@5qnHl+Q|AaFj+4sMz&weeR*)KElH%c>IK#hL(<llZ%^&mrqno
zTtZR`B(1EXs-~`?sbyejWNcz;W^QNi;OOM+;_Ca>&p+T@U{F-_`<U3c5Ag}<8JStx
zIk|cHW#tu>Rn;}Mb**jf9i3g>J-wr2;}erp(=(9emDRN`>l>R}+lNQrj!#a{&Mz*%
zuL}l%`*~S^U)c}q!i26179JiB9_jnKU|^k~0}c}&ft(c)>!~8rE1QQDY~IM&&mupT
zw%ns+S31Daw;e&jrQ%qAeE5B7KdtOvTbR!uTG`(h_K$Tz05mul=*5G>1cZUJOPX{q
zl>Z(7AMHW*=Vz=)NZ7zOeZWCyV$qHe_9&-JIKBlf=a=`~SPXAw_Ru!*yg%IQlrJ*B
z16n!r?||W*V(?<3@GYwdF)_{{_}8CFUQ=q7+}y9*xVb3516HwNMn+z2{P5WiKl$Io
z{9?cUcW?ezZvLOwDI6J<v3+}sNA;f%*DJB5dHU5{3^tJ>d$Oagx@9wq2U~+W4m|b-
zvz8j~-ySJ5P7J(#C_@fdIf{Uqmu%hCCdmD;6-5w!Z*blL6V=0Gc_p_+pZ&dfn-F;A
z{dGwkcwS}BMu7LCCXXN=6=Ztl_%8fw?*K$_5A(3_E~fPzkRDs-aFwwo^57;)_^Zqv
z04BUPy`lmy#NFO<pN}6az>^w=Uya@Yi=B4>flv?WS=JrU6ne`khz)}@LwS4LbavR(
znXf(c$^q<RgmzTf5&x9W%H~tqzOFnH<!I@(r5%Swx@Pc9p(S%APIX<azbIk8i&cwg
zn7Unlu>=OsS6w$7H>VuuXw9*TA_<vbi%{LkNsGLa{^!`*?8pPnYDWCcI@64ub8G8f
zf^oZ%cYv$j)%oJ8<lj%;FuLpPKO}KQby=WYST`S{KxipjobAIMsXCHrq||&Q^3s`#
zD5zF|FY}Qu?J&o@`7*qs#*pRXlO0@T%9<!vb}OSfie`a5#ty)V8`b(XB^zP%x>RXN
zER`pyd2nvYHM)7<>?N*~koC;I`y!nvQOb)iRXp1{NIvYN4JQsJL8A+uGbb|I_0|tG
zpE?nDrcb@#RjxFelOb2`O?40Q^j%)8fH-$>+{h?ajfZ78VBhDw7t}c|Awi>!Kgebi
z#WP3)H}a32<rvp{vZtuiTZSdd<pK{5wH3meDi+Wy`kWii!KKD!X8E{AAA3QFM4{~?
zt5*>FaJl-+JK*r<(N*Ym2OuwwtpM|v4;-TXk=8*kEJPQcG^8kOEID9XmW$UQj-ObY
z_>KA7YGL#8<Upx(arp8o)-X8_1Y5@`NO6<R<Ii5QjcL8Id3=Zal*B_#M@tjWn5<>%
zgKR@^>rc5)N*;FP%rT`Jpz3XxdGwiv&!5>{Q(j+9SPHjc_1WDNt>DL=DF7c<;m){j
zd-+qrcYuz$X&pkCpH+{!UzVp)ifA23zJOp9d}4XbnUR}bYj$%3Uc#rdl26E)AG^Ib
zbdGcjTJ&$)ePsEQoYQT%!`%T?U)%IfQLFBN5g8pE{nyxv;=peWBlEh(!9)f#MBlva
zr&Eu2`+4LI)};c3pN~IyE+A{-V36<Zrcho*?(}#^JkVtjDTb1Sg+0mP>wuN2+boX*
zDz7Dzz#^TZ|5286;(;M~z>+Ck!nql?<1J#9Kh-^KC-rwh889Ms?*V!P^$BK<qdd!x
z<^(}?Sh^ES3MpOvk@9Qk?9UE9M=xi_aO4(KfOrjKkIdqDsUoNxfhZ7_0q=Tc>W#Js
zI$$S+b_Yc1rGcj@7ViLYM*E#J$<f;I-%xQ>adL8%s5eMd2_*<OmF+1nY>gUOBut4D
z8Aimhk$v_DrxYrqCQl)Ri7tpFv=MwPXGp%P8N>F@MNcZiL8UTgF=`rWGq6%XGg?gx
zT$MS-0}B(eifYp$t^NJ|WkyTOwzuD?Udl7rRUoKQEgkR?82c9Y<@bS2+>A~cC0fP9
zRuJL_p0oL@sW<mVxOv$kcK}*d!}f0cMpc!guU}S!)Qp1Y?Nv3gzq|lgYaSpAO3=w9
zQ_D77NkS^FK)0j$(HEZzT=Fwx&?CPApiWoAN#%3@OALIi&FC!<anYQez$5)lTx$k$
z-f(NU4bqO;ER$^#&f4ljhfZ$`i_G{}&|K(B^@_dU5b3QvZM6YOV+UZbcKIaQ%!0V6
zD$WSKLMbrApWj{ughUR_-U3RrUSD#SKupG87cr1cQoZY}9x^BH%3eV!zmG&1AU%@q
z#n!-|{=_@WBSLm#fU4}p3poYa1Wmm1DkM31xKDuU#vKs8*)(UpPUq?k4f`&c&j>Gq
zfRrGbG}ym-l${J1$c^kiPk5Uod=KmBh*X?-M_BEUpK*D#f(PL?Z^TQu-gkcXQu){#
zjWKKLG$oIBoYQ!RJv^4#kT{^}s7RF)XFz#>yI)6=x=;S=W=J;&{~hpfkm;mmAU3%r
zTsCAm>At`6jc|-3OB+rDEgEZ~4OweC8XySdj41Qgnw9Oopv~Q>-^ZT11Jo8p<+7VD
zt1bzA2Eon3M{?It@AR_X`UbmlM4f+$^n%)x;12ktQ+*qxwV4Z~2hsfQ`1vpnKYO=f
zB=)7Di^rEpH$IGaz_HSlT$8;A{f7gS>!=&J+ZNnA;F8uNN3TxkJ<yZ!EzP++bX&rL
zsS6aluxAuGg_>mfRU2E$jO{%{mMlGpHp*^+Va(cUb+n{nQ@mp`-o09Xu122MQQ&o2
z3d}4yrGTLWj*O(QmaqD_L+!C`lzFwnGfUZix4oj(X98VhF$Oi!m7t|C*C`!wz(#A;
z!lX+}cCW<0#H#T4e4o;6bY)&F#WWm2REiMJ42jImICT}(QdGTcg2TcZrKqPqpf^|{
z)>Q5B<2S@*?#&7*EIIoJ;*Kf`@LUm%jY$uBXk(>TRm?Ed$zG;*W~T=MbV&bUX`*ZV
zvgW9^u4yar25R5oK-GpC>3L?CNm^bWLA`wwTKPMk6lKG!;&Jv<wuk72qV`dwueQF*
zvD!KYvU?7#4v*H>`V7wRGZ{&Jd>ZS-1E@clx-1F7>&U|AxFO_mo5YK?nUmnhfSPet
zkmd5#p6eG*jE9k5^5zjg1f1DiHJycDoEDd=Dn85iQ(yt`0JW)5!Ah{t9l%Z~U3@8d
z2>x>Mz9N~VOROBXiUoPBfTVC^BHT_rH@(!d+VOQ%P#H8wk7!M}1GZ;3&BlKgLlL4o
zzze*Vd~R{u;<yNvj`U)u9|ZEh_m{HYoeRgd#T`&!{W0ynNlWUH%rV9ly!IWinN^(Z
zxp}G1`p+^+1501;ZDGRMqO+id+f|HUmmnNiY+cDwl|pk0ipYY2Co;KISLiMF!$VYJ
zdAd&*njQSgfH!3UC?jD5@m!i@A%1#~tra!i1a;+xb|h>ywQ%BNmiyV2I`T@|X4!sK
zb)A$^&v5GUgfw1+x$tB=Ysm!5y0sG!&d=#+YcybKH}+?1Vj4Y@ek4$Z5_PYq&49DK
z$g~Qw7p46|V1&cAaCivLG72C2p!bc0>k0e=DT>?~XC+jed6KL=!TvZ4Cp-Hw=6KuQ
z`ez_NIAazrFC<tm-bH+m45_A|(@|O5M{XZ9zd{ZwLy!$&VpisyXS?RQ`tL9IU~Xqa
zuWno?pO5+3sjX>1UuoFAUWt&{Emo(e?TPV+{@u;1#1<J%2O+MRN2q%Ssx_riiu7iC
z6v!a-8d`F=m~K+0UB@lii3Xfk1-kQtXK|+^ixgcR)|NCH`x6p?dmiPN3y@M~7&`EJ
zOWiWp^i;h#kzi$kS30cYJ&|XnRV1c?xR^pd)U!6j#(JUr!|JA_fzNM)Z#r-Fy{P{?
zd-Jb<=f7(_svw&B{*ZLxng(80wWKiWa=Lx3gcpsUpxys`i+ukiKGPrjou!j7lOXGS
zW8DOi`2xqz`c?ww-l*HhJD@D@4&dE;lQUR+{47t{h)HAn`rSjt4>|+c&W4P8TNlO+
zuii0sO`V5Vf7{ZoqxLU4jwjhSBl)HrcVjs>S9+?qA8tPOX>5<9o};qNo7|6&E`0nl
z)wNt@J?e{p@S&k?uhK^D{vCj@35{<IP5&O`{K`NB8`nR^M0THHpaF%`a>*6KEEp=}
zLB9d~<nwjr;f7|*cre&XVMgNe9M;ZHQ{)aJ%CL}z0=ntBI|P}oOwl%Z@`Mez+a!u4
zY<ZEjwY!L*Q^XPlM60Y7;)mjqloJl2LY*MQHSL1;?>?Qtyz9=-P%AxT%p7)hYOj)T
z?u~y^9z`D;3Oeoql8CS{h9~<vh~7R5#S8z6HmTc`Ilry-cpxk~mnp`E2lX`5*F|-R
z4CwQbkg$tP3NwY{Z%y#K>C7scp3_T~8|^d<B3uu|OzvHm8f?aX&$9+|_Q?e&r(oyQ
zl~-f`)i`nHv!|6~8#Qp-1`=?GxeeTzdZw|ko$EXDv5?0qRbhiYMeQ$5c<1T&hlmar
zc?Vt#e|K-dx;s8JqHGd@dNc7Rs5jf1+7HjNBl(Sp%nLVyt(1&~%AqB{F5!htmQL13
ziskfI_(-H&JnIcREp^5v)+-uc+}z&jD=fH4Wr?DZ@a_!KX4zYym@C?T?#fjZ#|{Z{
zO(`(S(5T%1$m5Sd?({$-mKsI`p8>=Jz{;|NCM#}iC-ylm`tykVd2upK7r2R4BfX3S
z!jOr@UpasA`!vi+>wSgppC3K{3I5%17X^|~?_g%iyqhXB#txsGbE3|~=n~j?+sU!Y
z(JntKj-5e_RSAhQeShnS#5wJLs?yYPBgJX!;N_lvO;cQNRk#IfWX&-WTx3u4n%#hE
z{h?!|QGQ>x)Kr;HTls*DhtSHKM1OULcPfT%Q%EZodRUMh+_*D7-DDev%e|Va)T^iO
z-^C>M23x;(ybm&n{)S~UjJm*@Dpp!I$;@$MyFb;DO?6bbXZ50Dj1Gg=JuHg08$QCB
zI#PrrMqCn3dh(EC8%xQn&NqUPD10z84fgLUo|8*@E1Z4@*ff1eySEAToiu-HM8+$k
z)+4CmQzZd3Fcuivhn$~dC^TCh>TcWK0r?dHu|m7xw{DXTgH7O)Il%yb5lq@yp&YUs
zi0fO@cJNT6#byc#)NvXhn^|S8Tja@YQ3-8~L8{giCf!b<-tE4NytJFs9l&%4M8@Zw
z8?Rg*a>`@%R;yd|s}<a4FK!CIySf8}D#u5}#n#tPINoC8dW&PL;Qn2qzh^Lifl_xK
z>J%k!_!lpK9;Y1pc!6(CWu)pa?qE90*meJr?^M?7N$TMB;i7AL59UBsLMa&cgYABd
z!H>?mcJ)@Q!z*4qZbepKb*uQ~s#l;!M>MO?uP(Otj9s3@MWN6!d%1=8HLQ(pj;`SH
zp7t|Y<43-#BQFn53Ec?>ZVIcBe5?B)GBkti>%2U>^wXw4#e;t|_N!1)S7(`l^x3^x
zCJ#idCuv5De$nZ7h$=O8qL|Y{+|JG-q!}ZE5%z{f-pVCA+BH<X(%#{{qBuDMPQZ$#
z{*JUyC%O=gVbeI%<mlqiQ0&_Res-817{xFcd>(x4q5C_J&YU3GG5J|dv&UPM=!|7$
zT@1ENYcC5y%4gpA0wubpaGpU+PY%DyugKKXyWQ~hjt(AZ*49-fC{wP%#Xoo(5T)vb
z1;X!mG0^}jzNy!_?eu9Wvgc6dn+ssJU1y~gqZ4^ft}A6F|1=2L-ON9I6lu=qP;`uc
zRICq*7Igds@9AHt0HYv$4Hmy$1+yxTW849xBjB3Feg*jB-^pixB`v=*2Pd-F0z8BI
zK%=3DY7mA8K$$8Q?r;swTHDIfm~A{Q7TOk@esCFQ<>dyZ%WWB$L6=J+aNy!l@2GRv
zVmE&@$k<vc&&(e`Ta4|rP0wk6qu#OqWa4dtXLwRpP!fn`*HAE^`VJ`0J~;CwA2XEM
zlYZ)wdjcOm-?r;}c?ZyK@5q^sopXBJ0e#o~Mvv6gT_vk8ph|&LS%K@gWo?Nij)F!<
z)(o4K9@J&0Quv!I`U~{H@Kcw7F>eNlQc!jL0T@oc@J6T7p+Yph!{G7YVU2G0W+IE7
zPrjU8*g^%&des7D4WC45zPx4TgSr<4e#vDJnxmPucmaRjKHj6&lvGw8oE_D5;pyro
zXbkx_@d~3k#p)~f-W^~bZMA!Iv!KhKy0VjSy#14nqWuZIOa3|Ngy}zLzXf?hW0;r0
zU})+zRx(d2y)N(@4hMJ&Hp0!wnbVX_rE`XGL0nq!R%)8i-8sXbPCZnLQ{YA<x{#ur
z^^stclXHGj-%FdvQ5d9TrU?ZTrJ^y~A%c?y2yZXlKXrj~bc*JiQd1eb)g;s?BT(wm
zyxzQdV;&omX))<^=32Oyowz5l|Ih)eyyDHQgx|ByX@pc>UNVdihaC<a3P^C^rA93&
z&OMS|9YK2!+pS`ggl8+xsD9sIb9G9RYvDh2x)P~2s<*jrex4F8&H1FQV}mVs_PA}9
zCqzsK)`yhlk&k-h_<B7b4Q!)dl9pq#g#n^_foT->kWM~8npW@~`-HmMzkyxye^&sJ
zhfn|mfI7&3cfhDC#w?9(#PS*hPS}#j`IK=adweBWvFNefTk>bo)L05v>w#XX$!klC
z%Z21jE&Z#r66QYQGNz69mk5P2*7p~jgZ3OS9&~=v4sqdoZ4z=n<<-7NQyfJ_D1Wh;
zX;V?4LwfXE)h`#WVbQ7+0-BWCo8r62pB(%0_?9zIFfn$n>%VSZ5{Cttt!au^3-nOh
z9+9dEZeQDEDi+$Lr_;3NhkQ@473tz$l#YP4BT*dPM<?z8ZRmy5O!!?m_y;KX4q6=t
zdiW2fZhZEdph>g%i?6><iT)Xl$cS=&<s7AHAQ68mR@nULHE#0>j~WlL_i-LwN*>4Q
zB(sf0PUTIz7~bB*2!0?B6?+6C40TCv8V2*#YvF@%a4Yp4a4+ty@TcphvYV2<)jeU8
zt!iV(ud_${^IyXzc4EF<>TCpZS|iK^gOQ)rg-!$^OpwUPGe~s{dwenpNuG7Q(r<j`
zg;E4&oa0gg8JNeZQxcL#>CtLO+}B3J!xzMm&dJ}zg|wIU2u|kQ3w881$(6Oj7?)^%
zxNn|w%4f>wF{UxQ?OOjFR<LL?^Bb=2)7gz$#TMQ*V=vN=7H0%08F^-Lr}xY}UVY&x
zK<`Q=Wbh7JKFu5IwM7ZB+0j;FL7pLd<C?J6K%D-<EFnA31~;PaNz9v4Wcv{cR;<w3
zFRfUCHQ0DccGXX+Ij!$+CX8r0=#aEtN{#1p4%A7hxWYECc(GmohOYj0aL9hLb5YTR
zvCVyY7-(Eo?V^a>t&bf+-*roTEq2IMvMOUKtWbBQQ#pCitFL-=kmU^zvlZrHF((E8
z@`#?V?dd2-eIOG9N1u#YCY}qqpK9W^NI^)HEC<I#C;?4+nk4mV%p*#XF=Ogo?#ukA
zj%g88Pl>`Ftd_IUP&UaO58c!5$aao{8iJsA6|<U*Yd#295fzbNCL|BSa7u@*geph2
zT~%Yt^(VF7!kXgwLqJM+ME>wE^WY61JjIcoZB1azs$dZ{ekN+qjFVu74|&N#oK$MU
zv${@Qodwu_W<2=po#&h?G4^W2{?NnC(YnM+ZAM^`0t-LzwaJDZ!A6aN!g+lgx<THK
zE=!{PaHWj%5(omhowCo5$Qk(3lSE5K<en43r|w>7l(79Z{~dTt7onOhmip{fFajJ8
zQIs>$VTw<A$u()Zxzni50;G?>yLVLVOzDD21wa+Af3?Y1Vh&9E7}Zn<*~2ruBQkus
z6f*w)FSI7MDg5JB$9qh;)iECf=K_tCeR!5ks=!$s#k)a4^_9C04*ll3)Uj@zqYeBh
z3yWb??^EYyuS%h|;zo7w+@(QurBPs3aJLE~T}85btNmQKE@y0Uu}S)?{XA2jRJpbc
z9xh(or)~K>!2!!Y0x)jLPmXbjnzsJw)?!o6MSLBzKgC#%#nlWH|4@zoIxllA)GhoL
zoNB(19@MLB>tr11<5Tf0h$fHbzl()`?SK9S{QNhdbCnZWxiN_`>T&B#)>)DIWS>~U
zHvn7|hz#J179r8JF4A&#K#as_Z>6lZW~7wGlP>-84!JP@DL#To&lc1+oUB1~1D@1o
ze~!!0tNa)&7a%jPKVGp#ULzHgoF0Rng6Hr4FhXhwufuiV=d=IX<^MYw+224|WRv=;
zN>jC!6`BMS%yz0b*^({%>Vk4GORkFO^Of!8;t1_ZJ#f46m~D8oz_Qooo2M<i7I0l^
zM~mkU6@A*YOw7zi8b$M}RkNVvR>xOo&kWuOb^{LC=eQ6L{%6t88!BSycQog%OdR;F
ziNjX;O0Q<WMc<6Y?1>Ny{#QBgFK}xb)g&}ymL>`Jv=lj%7(Df2ow(e*=fEkeyw#qR
z4U<->atEypAvv!veZ%KPZ@-=;UVnou%l7G;%gkgAOj)zbyf}QE9ijw;JAj!VKt%g$
z&yGn>l1-ND8+ErhnIF~-P$adzOlrfJx*E83DPqKEb4*~m@D&_5ZWwUlX4<&KFrO)~
zQy5KMywwSr;1)FZSlpNs6Ayis?x)xcyhPEWDU|j3S77hI0KwA~q8j-!Wv^c9v#TY$
zA<AKg*U(&qok!NYW{CA+XQ47T?LMsl{YvhM%L(V90~F%-;di<N0!2F8O^j9tcH#-@
zwU#QN;VWN9g6<t)?F^P>GjVudWla^kA#5|<-^8-%Of=>8b20V@Z=3qhg_jpyzQSv9
zo1~}6(C8Ohb>h1jIZi`J&);@^%cmCSko&^w$THnb&$OpCBuPEhS#o7g_wgaN>vM=^
zCHwx7p{CQCmmJ<p120p!(VvppVz-Wjn3GIb)v#o?*?DuPK2s(Xq$72M?tYl1<kATC
zPI*@o|B@4oJxdF*ajV<eCrO7K>aqkGJw$;Kcr>Ur58BAkLTp5cALDd7`+Hc3HB{93
z>bk`je_i!Qr-5{0NzHT~AQK+t)dpL2^>wf3p8Kjy9)HE$=H`w|_!1Yeq@l&TYiVg1
zh&b0hYAq2BoGH9v7a{Q#rzorQSHhR-COZ&k1aP3RQrEHGKKC~q!G8$%zgVEekyn~v
zI1Xla@L&ld6B%R!9Rm(&W_n^T@_zzAW-8~8chF>ZCc@K*2VTmr7T~W+*=uBuJKO_d
z&gSFH?ynt_ZY*St39vF$2NyW|7&VSdvUY_-EiT=e&>e7FH=zu&lNUcP$Cq7<Le*;l
z>x4NwTPq%hWiaKqQJwO+2b{5X3%hP@iY!a8U!3jedYhP7=RmqF%uqzp3{c*&(7Dpe
zWyB^jbx3ec)PKGMrbN<%=-FET73uyy=P!_OC@_o^y?IL#HS|0PU1nLe=mug4auKMq
z1@R6oIvP{kY`X|d6{?4dXEJt%M7Avt$jRe<M%fAMw`8*02=<z0J)7M++AA?igN33S
zGVBNnP}~_%(=EtaT<jy9**Syg=PySs;1Thj2pPTxGpGsWEL&QgNu)G8r(2W=t}Za|
z>C+g%$A&+6*Mr_Ha18In0h6^ngspL&%&gnW?-w=-rmE=12EA8TmjltMj{CfE+0thI
zPEDv8LeQ|7bcJf$qN;`76YX|GmYD{MXJqfJo((RN!CCQMcyD~#<n&<rra*TJ+5ebB
zM~=X#gov^4-aUOQrhsVj&$?7>@Xx5c0z7wYgx-w2i5*_lJt)*gx}wA)ZgWd`uJh5^
zM2jG+B2IbG>)qR3<ZuBj7DTywL{rv3*MSY4eZ~0Zv+6qREt&S1u(l_c$iPhE@U}UG
z-fAx;B$0CSwVeD<Q1P=HEDb@qx;kBgyRGANGFVI_DC7p842af}HYTh0t%`PwGHdHn
zXvjjWT#Vl2XN*TTXKU)sY^&utymCUAsK1$c8NtUr97Y58K^}Ey#5wz`wvKD7T6~pQ
z9DzZ-c1zNQRO?Gb*Vx%87eN5`WD#Qr3Uf~<yYhQb917&Ulp!<5jA_g7)q@i_Qix1w
zq7J+k#?9vF;hbk>?)M=#pRx&bN(nCK?de`6-cj?y#|FXEqFGe&jMIgS@L*p;W$=%J
zW#taZwq4IEzkM}ov7!u=&=uXxXyobWnL$xfR5yMbMa5bXRhMbb*53tSC|L}hoYfs}
z7sUo*ToPCKK`!-GiM=SzLkv%LZsp#)%H^!o;%?}on`QgfF>(?j+@R%Lf|+h!vvgQP
z^~{1uiY3<<z!0Igc8p7ELbwG-M%m@=SLg-b0n9_gUm6}R)4Frw^vAJvo;_T`BqbLa
zvO6g*3!<k5^~@(~GGcvF6vrX^O_1Ps;E%roL>rjooFzuq6#7W!PgoFM9+(`X#6`1<
z8YKZU6J7wbP?1lVSrv;Iwx{o`=TqjTa*GrM9Nhuv3w=ldg=&2(x}oS9&r0AQak0|3
zzBME#ag-Rg@r_~)Ix%mp8WvzSvB0Syu!=1e^L~Ig@(y@u>mnR#HT)Cd<Bt}q2q$M=
zy;(3UUE$(iksWMsQ1p%7^VQW`A~tA``%GCLpjA6Ao;p1@9aq&{C%fpWJcx&0bIrT4
ze&YHZCJaiEtblH&tuEhD^||$)8yL!NtkS=dy8~*&ndF%rox9?{oO}%mXp7@gyJb04
zkL1AAlyR~oc;%BZkJsc}FkJ;}c;8OPVW6mSuOfCzL0kX}y#zZKU&Hl*Gp8;wp-gq9
z&m=n6AibLl{HOK^x=cBCe+5LCwIe}OxWLJb=vG_ZojT*1a>61d?zrbsJn_a(hJlfI
z=*tH@YawBhTE~2YvTzdN%9hfC6fM@%vD`AReU#*toA7TeBz~@H{QDK6LQByO0=0Mt
z!OFRHTK#u2{z}}wGam9|mjvUrbw0f_`%K0^7eHYf)~}v79p09n_JNElM_=iD8e)4{
zy07TWLrs)=TWKbj%2gsh7pH{j7+X@t2d%c_UL0iw1T1l1h6LfQVG;pn5vP|9dA}~$
z^M}U@Lksp7D-|5F1XIt(49Y7*C4g^4X!wNn1S2_gflL~f=jEKT)dPay$2U&Ul4u#%
zAL1_mT=C3N*87=ca{EI~sCMxOP~o4MX-X3u9_-px+Lr3_RgK;+q)O++s1%wplSs~~
zzZFy{AO(Y4am?-j<c(9(Yx+CjJl7R$dcks&*gv1j?fLb@<GQg?NT+g-b7gYUwW!>J
zKx3rRx5!+}p|K&(_z^4TN*}As_3F`En2c9LBG_ml>J@)V5ddS9tiWbCs!8M26~Tw9
zwOp5;1ttzGN}W8$>MxPp@6QE|#0YdJ!?YMtG@LarqCGxvY`Baz&$A~$4Cs4+-6@^4
zVRO3&xH3W`dVCVH$-zJ~s}Zfr`T6-)Ey!%Zy>+SP6pam2ltNI;l0C-&l1`jz?Cs!b
z->z2Z+hF2$ROW;aMreg4-Z5lnz#3+dZ%2aV-khuIX)<z6-57K1u6EqL4W_Qp_Ya%-
z4XI2*!qW3h4=fbSFT<*Z!jxpywvlD?YQsG$8%!lXtYH;8VL9drU=C>kI^g!Y+W-#9
zD^&TK>Y9eI(dpvs9R|NyVE2fHg{|Fc32s2rZB_e2gPz0km9`Fbzi)3#<AmoB)oiDy
z)OGfJcuAIC;}|n1Dbd!Ha$6|#^J*7;{kF3>Huq%q)JSey$?J;MEaUUmF~!M32OZqm
z);`wy;ND=xSW_`+X)h|PzcCv^2!K7B+@;X;>LBij?oDsaAvH3GJR0C1g_^wWzw;Yt
zU|oNYZ4b4TfuG4IId0S-oxNZ7Ec+>x4wG<SNA~P^hfmhI-2qs0+mFJt&W$;1YjvMA
zSf6k~NfU>p)*Azwqt}!<GX(?g8A_i&tXm|;plMveorzoH<$IpU*I!Qes9Hi&7sqZ<
zSp-<Kbe)S^r&JG#$Kw=ouH7IpkkByh8dB%fZ!H~8z4<{nA}ELe9Do9SDgeYlR^f)?
zzy%JT1V0p)Q_%PWIOI>b<<F|d_j7P)3YU6-aUI=+4xVP3{fX8E#R-0(v;9+9i>^3}
z__zuyJJAlIai%Md;UHeo1;5u{$+>%*ZUP%-rOI`Vy4CTPr0G0zxxM0w03y4lRWDO+
zHHk_o#oOs=`T2ZLcS6@Df17WiVf5S_T03qzbyAK9`LrF**iX$n7l}9il~|Gz39*sH
z3&rR3a!j@&YX5Q@RonFNv)LMQaI|-<z<nvwP{j$W%IU`;W_+rlYeu$vcAk@8<kMiO
zho&gX9OwmQ+O=>Nl2?#lADd=6KUan<dZfi_!mVOp!Ssf|Soh+#pSO^n%y8i-k{i6_
z8Jx1voK(mej-ke}`HA5=#klMtad7UdHNwORKb}#lK_yd4KnX|Y-gyNmwY?O8X$zUG
zDnYDgRnFYz^C*0U@u6Bxp4VVkO{HrMS`==fF+4X}UehhgoQ0Zx&)DkNEm*G}8lE*5
z<I;a(2Lew|Ywm!0=9??)gN!}s`5pDHP8XMzM?@ZY6btY7K@Tu{O>^}EJRf&^Q+=-T
zoh@q-epLo#%R>2uC_+=DXGT!AtmQ4vceX5K^>(-EY_F*`yt3<N2K;pv>?be&FP(xm
zx6>jkQq_kJGZwFh3hbx(Tm2}7g4XSmmnwNC*8eC>UY6css^Zo1Qy6g~xDaUhXP9Ma
zHqP&Eq0vr!I5~x@JuTKFs16l6OU^cxtOo4|-1g>NZ&W4=w{btJ4$MD~RjsSn)vwh1
zu0>Dhel9{zCXM7AY(wA2T|fyUIUVVHi&NJ@lLQ{w@usm=3REoA;)<84%3sw?)60J1
z<w$o!5hu|<<b(dp3jBFr1U1SeqBnjk!spL~mr0<77%U8aLZSK7O@aRXq2!(Qvq8@@
z-aXyqOkKSNh(JPEEX2O0_nnq~ob8hPj4odk>z=M=`pjlroeKMcGo|+AhI{M;jws^>
z^;t%Lj0WOP2*BQjvuh7Pxu1lx(i272MK9~B!Hw950TP$K9UAFg^cq`}d^(@1#&-Ai
z$dD-LontwXHUS*|pLj)o#>@T}fb`2(TOqc?lQ`csp9L8grs8a!J}eeOb)Fuf>(#J{
zrUDBtzQ^&kpFc=6;z%-Vj7ZwO;L2KN6zz&LP;FEhRW&%;i^*7|eX+sZ-Q|vP<d`EX
z?CMI%s{%U8a=NiH+!hj^U@o?LJepnO->Q-E?Ro5lg+S87n4xtr%CY?YO-BOR4+%4&
z?Pn+?+O@ncpbcKZ#=J|6fFcfb0zoc{v2-8&^%j#Js*PEyL#-Q{JXOv<1fz?~(5HF{
zT1`9F+soqEYL9=nwDHST{@>pJ0L=V>G6*zR$Dk~&ljVZ<8S5mW36=t(nSpO&U@DaB
ziFgCv4guls9eGM#lyKa3$K5^v)8Aq_-vP_0S3-6-NsDGz)MwoLo%0P=GfAEfTlRmE
zjsG=#n~~yNl;)BXNfC&sIkyC}Zwf2Da;?2OeJyNX<!B@7^hV_(K4OJ(mdxSO(Ecos
zYJyW9S8DPS-;I%hxNMuL#!OI5^b$%vxQtD~7`C)j!q}aro{ch0hO2)#9%bQxH#V95
zGt?ZEn7c<@|6S{Ejlur*;*DXm?OKJMR|++uh9TWTQZQ@qy1qL)Mr?Tc)qdpuoMBWj
zir6n{i2a1hlji1UGGF1DJxC9tL-|v7;J@s__t(G1CU~acM#<D(QwWGU+ytXH*l~u>
zC_Bx}&oKx&+Et&uX07Ib!>;SCt7x&xD5YuP+CJqdv3Fcxk*J}DpcW<(>ww#Sv++86
z?qM7IOSGv1XaM@VwDBq0Rbfo?F=wHj?h#H#Huc)7YGx-YP$N~0qL*o=Z3EH6PHsa5
zq_m%-0GpHtRXV8ceEvR{R&kcRl|kQz#$K65sW5B>s8Q|1li1o~yjFi&upO7+qEa<n
zi7u{lkJ<cy<~v0_=Lpe@*7=7};h&lBzk%HUjR=ng^zdM-Cpcee)s1Ib(u<nyuWEeI
zeViZd-#zP3+=+qTWuPvXm{|5iImC(HUFFc}-La)opX><L6_Na-WX2hBvLlJb4-Z&F
zQm%vQm+_bdT(Vm7Ev7<kj8I){v^!I@<oa47OQ(`<B`ehr!2U)}pT{3{)QMFq$Am^%
z#A3B8gkhVUdlOy0wEwy}$%R{I`Kh|xoBS;aCL%QD*l{HBseHdGYun~|2dtTNZ(PgP
zgHuoRomZtYG*-pM6#Z2309P2f^q`GDmiPQp#(sSN=Vi{wIa)PQ0t#MK`8}y+6NhK*
z?ojDr2qv(iOq*L3nkax4Pni0&#xg&RTTZ%ZO2>r;0x>V!B8TAN&J~bSwk{<KiY*|4
zaZ(2h@r7n*&-g=)7lvgiV?(Kwv1*p9YWi4tb}_8wiWS^B4}4-NB}6p7yi1Uim-95=
z-URZ`Cb{a8LnhTfQJs#_f+%ewS_`4MN`hp-@c>%?<O9AQp+zoDW0%Y~P4<s$HN2y}
zlh6`Ie$9XVgKX~nsf5sEzYHCn7hTvWY&wbIZBqw@^KR{GuF}Ce;Z7br?P_n1zNYFu
z{i^G++=OgCdAX1P$&MP`qUtrMwn%Or-k_wDCtFRU7L?DjM3J5wbJR=y_&#6fs7?Kx
zS?`E1gR5wjy=uA-8pimt;w7ws%Wo!xnufCMq1t1Tsacd2f{I2z?La-AvawL)1vKMf
zs|;pM(T&{4zR+WzID;Cca*n^uGW}9&{=CHe(#!oxk46`plRsm6YU$*`*YRL0rO*S|
zH_MzYduT>qRvlRK)xHURanWqc>WoFqP`t+DolQ^ArODVf-_*nO^Nn^-mYS4V=`pg>
zkS=`o=T`}68+Mq(YX=#`fr8G{7L4+p=a%B6^261LGS(NCHog(sR~es$FFWhdh1*=Y
zqTSD`2<&msFtA!jDVAb%eXr$}y;FRu&S;Aq;+8_J>FQ0J9^Hy*kp?hwD11Y7OWblv
z*tgOdIwn~htbzuT?~{^e3iZ>AzPk|V|Hy^<%UlXH1r$NRNz8K&5F0U+I?*wZFO7Oy
zi@W`Wf<aJ3pN+20yv{nXh`M^%UB=2d+|1>8=40{-{9#1PxCy+ipp6Kj%Fp%87pcC~
z6HJxFF8h_i$&zx_b-1yiHrO;BtE7#uJA)$|_DJ_hwPq$xt^|n+omjT-kCKx?w!Y25
zoUmrVoDU-~{{YYYKIbFy0|h&Ym*Kb*TOrqI?@^(25Gb;axweFGR6u?u)I_~(c%O4~
zImjf1#u}k`%rT>O6~^Zn4sVE>GmFp4^;k?UJ-oyKTG;Xn&fZ2QAEroE`W&&s^Y(f0
ziGK4N;X%6f_~k&VMjh?;hefC9yrHaDiC{xLz3<IjppVY<5(g)v^O};(5HSf#`fj-|
zJk(z{`R@DF(dDUxctF5lB^++*0pWTmnUe)tVtM%TU*~84eox7NcZ?r6hr&NJs}v;v
zkd>?aF6F2zh!>fefR%IsKN|$|K)NcM=dmHVUSGFb7kn7I?<^-LHcb5T9sx6kvuB@l
zkSI-cxJ7RpqPjqR<aR&k+Cu!-{ha?|^a@4i9iiy_LtPnvC_0b$6FP6a{h8Xv8KUL~
zG>u<!VeF@!Io`g%A#A&6V(%a*X9q42)|_)My&~Hc9>OWhPOr`Pn1BAaE}*9Vk0g<K
z=lQC2Q*ByiY5T#?C0eXSX>U{env`*pD=<5V)D9`0u@dBaJ>1$-9#P~A)3V=H`MD2B
zQT%DzAB95xo!|X=*yNWX^}or0P&*mI)Z8=;e6qFwWra{fVv9;xu5BSne{DDL{d?6l
zQ=J8EjoFGb<t{X7Hl116T8YmS?)sgKEIm1<wQUQZ+RMXMSf|~feM<#!9iEtn;lcLP
z8{CeQQ%{~}CGy@wbLPG9;&op=xm0&sE)w#fZX+;GwpSUlv?Nl(=0%7R;Sq|du9r0{
zxPEjIuchiV;>mm+SW!gLCH}odhCwk*<HSU$&@G`+94FYazT7ZM?Ws7n65hYFyMJcx
zf8@$Ywb@48XJZv1v0fI&)YWySZ@ZKXc=b@HD_fiDObKINo0QQmRqgcFeM;x!x%sS=
zQ(dfN!#K2QS=))5!6V>IveLj*7w=Fqph8BFu1kA(a$)jusEp~6)oU~>eI`Uu7d1}C
zDKnI~HS8t|#lHlauA{Z;5^TG~FXHJu&p1xNk>Bs#Fx_ARg+%tCyjyrEBzYuDs(+Yz
zgH$1W{jd+pwd0J;`@c`(_~|APKOV?<tD1y;-1@{#Z!h$|@a|4l)%;B~?KgHb+x;#G
z&mEu{VmBtZ30Ik3J5=gCbn6Qty}kDZ8VM%|b*u_sVwi%@spFvC8>Vm>|B)F3lluZD
z85-gaH>Dp3K5s9~TAw~I2<I{kcJv%}wI1YOB~n1L5VuSH)pl@r1b2duNrvx$19@g>
zyTmY(<p10h8!gMOkQ2fi-Ki@B(eNw2Ez4@KVXu3C8&Q4fX_)Y&0nzc0>*&>^5NJb}
z;Kw^)DMjxVly*NXXyb<+{f|vqsjLsj-<-tM0lzl2v#zVIn!SmpTOEj}m`|-BN&FUe
z-i4F(gGByxbC4f)b>x3=&ZJeM&8ncp9iU!$J3YL5J6NY}xj4g?%FfcowCj-(e)6G9
zkbS|6UyDBJhi%K))n&f}(oO^J03R3d#S<?+zW>;q!g;IJPFK#i>Uy&TM6>%vdr(l;
zu(#%$u~fCm9#Z!8>Sxvtz_CJ|_lKR6{81?Zd;9<CNjTYa&O?e#6MNt)4yf9f7st*z
zoj%X$?`Q3+CpBHYF=iclepbZ%Ya0M;Zv8|cul@S#&qWPZ*^qA+hM{Eo-0ucW<o*v)
zjB4WY4a24<(HKWG%ZC}Q4(oNl*(QBI=-mmDTqk3GA1e5bLOM=rp1GP9m7Zx$(CR~w
zkLByE<oMEwx6)l{aMZ&nMR`wG%$QNlT*D#}iV?B4TS0n)w?^?2!YKY5%Hsn3%gK@|
zQ~R%rbtk!EJiDw8{DC(gOl&o(iQU3zVAq3J(D~b0=47&a*ozuE?}=57KI<U^V(2OI
z29fs`Dt1s?F2x|+{)~+KgF6>aFBXK)aw6g)plJMakyF^}pPNd)2Y=ZOW^{g$IO7kR
zS%%2gNSnRPi!&2KVvXzzf;E_&&;-ZP=?t2Bi*tlSjYKV0R@QX3>)2Zry>TzyOFv5|
zJ?J!sOCs_2Vt_bjNPQk3Of|D%7(+B%h#mb<&E&+KSZyqUL0GR}Mc;Le!)|uqtSnzT
zVZA=7cso>JR@U1eA(nd;EIBGz_|O6kcSvo(4dEC2tHB+6%V^P~)u%E-QN?&wkt!E6
zhvE1zm-+x3hO1x&ox1IeMciqFv17|Al(Fzj^4c^uD+Q-PGkwo#5xq?AS4)c8x6+i6
zPBgJk0!#M`;V?eJ*6eu-lCWokR?%JQ1(o0vyHk427*u_61{af}-V}7Mrb<s8h)>*X
zyxpT7PcZpV@sTQ-3&ikD0+XK!5tU5IvjS%LA3=cfc(h_ve{qbB1;XHRK%IBibEgZ{
z@r}bYRCxmFkA~t7@GNNT-MmN#WCTOEi=#vDz@Nk+$@1OYh(O^RFYl=7W`QvFl&76^
zZxx(7xtWDBKFWC$VoJVZ5KlCNaN{S!IkF$6u)lyNMkQOR5A5-*kNZd~`b=CNPqdn&
zQ$E*R`bI_Bi${rq!)Mvt4Pk@Xe*X&*#EHw>;j?<%Z58X_7NvbdFsK<H<aO`+_CVOc
zS*i`yhZefhlUmt<IHlI}nV0s(OKzW%p9*{ByLRSt%s(NNb?Q*7lgMDwwA%Y972_p{
zWBd*28GDb=hqU_v^5cx0+Q`<0CZN4+!YDBTlEoBT)>$0mT5y2(R7-6BaK$U5RD-gz
z_^AQgCv<Z@7mL0(`#c>E|7b~2u|p0K8_YF$&bz>d2;Ct{VL$cX+0DXPscS`R+j8@S
z?BK^yuwx7a?xjL3uI$g7?RHIJF0F&@b+k@i?=fY5n@8=8*6zVb$jE!P4d)$sW~c42
zkw6-FY)x?9I(jRp+oE62mL0E8Eno41SBb><`4@@d7Ipkija4tf`@d-GEF|m$<5sQ}
z8ReTsYUAEO<a4C=k?QY_KRibP;$9r$1}~cmP2@kfjMCSr({Awl1W@fNIpeKesSUg#
z#dILuUGyhNu4>HFK!3Ikky<%=-&)T7@?C!yPnVrec5>eI4lc+odqLF6e9bT6(fv$c
zUmw9u*!f=~8-Lk*b<#uJO~@9XrbXf=S+q)Zn4RHKVG;h|7ZUYByG`D-`>jFWbf`5a
zDhk#$5Uvv#*GnwT3?X%W8?Q~Ld)52fKg?yAxu4}Y?`YN6jtZ#0B=?Jxs~HOcjAY_s
zxyZvQOVX}4;YMjn_zt;B3MD_D=ZH@Fg+vJl6qhirmrII6Fe68NS=Vr1JaFk*(k}!R
z>zbWu<&GUEx)s#R-XC=C6?;A!YbG&C8LJA>puzuAI_bp^q3?Itv)WNRcs^SEF>9Oq
z!FDVR)+9b-FE_-1+)&ckNRan<Pb)5kZTwmw;aF>&D~_AMqwta1o^e<(O+18CLP^7i
z*fnI01QFV#xgGSa#*@q(SHaxQoXdY3RaQu5B88Z1eL(ahX`RrUwWNdQQ(o9i^b(AB
zpYe3;t^|_{zv`ez)%o!R+dgRDAbRD`9^KqUM&d1Qo5Fmk0>N2N%+aC?PF$wH7&74;
zZI9gRIa<E>5Vg(cJw=_y3tAO)UGD;EyP-OlSe=-R&F>ej`*S~i*!@JzqZp}|hj+70
zuo)>$@dd(Em;*Y!Y?IF=su`!YrIH{g$pXF>HJe0PVT#tiCuW!m7D%^Y?J0BENsWMT
z;T$U7MqQbkUir{tzm809Sj}fRGtB+?E^AmqMndUcf&@XUiz<hpDT?8bagtZ2Rovx+
zlSdc&Y_nLLIu{I5uS<EFNtt=_FR4QuWcn~L^vt;?BU*CjU(CF4c^V$g*1dy!PY@$k
z<*<)!m6!&uCcV_1xnRx$@*?4V9+tMO4fFIVEUHLRC4QNl%*$%LWk1XR{w@K3`5kwD
zxQmbGWpSi+K`{B7qmmh>)%x&0GNUmbB(lSyJ3zOBCJk-gTC94iPzmB<;+GLm?eC1o
z5Gj>I*DFOStRVTMl<M3q`-O)2X`DJFrt7r4`OFT)b`hdaR*bChh+l+<m`}eBLIW2M
zt%_0%%6O<Ua}rCOT1xa3lB}$Vv5CXMqGV&koCWs4M<Ln;Rop(RzEYhS=dP|gdxRUI
zg-`B)fcwN)saCm%$u7x<U+oJvOg{$Pbc4FEIkbUt-^!#a#Ud|Rp!0}Gz`{@D<5Q&c
zwY6{R%-8hU+fH?U11o(rP32*|ropcqDfwyf@X_N0Ha>OITLl5lzOCJi`$Z1u+R77_
zR+f;NrN&5<w^kT);#8|4QSyyXbDWtVZ5m}0k0%_ocM@p~eK0EfLvXn}C168tMuwE<
zQdq~J@HXG1w{>A-;rjaa;L*s3*?CBfCgmJU+d6_28JZ9BlB){+@`n#wB@VfZ^8-2a
zm2<X(o0@Z-Bg>!CEZb7B?_1ag0UP@bq}dC^vRP4UAg}|@awHYXuq*)u6kcD2f!Pq0
zHwXf`KZ$cy9q|8HxE<Fc5K$ZD6mIKcMf}W9B!dd=Fp7Ul6iSlwWE|Z$YAPSMpA}U7
zszmLVJ!bVXu%Ybf^sN&DXoE7Qi`W874|)5|HDmcflB0xneZUjYYSQa6<HUi*+IPnE
uDO^4vU7LhsR~lGUdLE7|R5yCVd)C%i9RwK5s(5fmP?Eu)l2W9*ng0ViDsbrl

literal 0
HcmV?d00001

diff --git a/validation/lj-mixture/reference-data.csv b/validation/lj-mixture/reference-data.csv
new file mode 100644
index 0000000..039dec6
--- /dev/null
+++ b/validation/lj-mixture/reference-data.csv
@@ -0,0 +1,24 @@
+t,x,density,u
+1.2183,0.5,0.8,-6.991912
+1.2183,0.5,0.7,-6.28399
+1.2183,0.5,0.65,-5.86756
+1.5096,0.25,0.8,-5.87326
+1.5096,0.25,0.7,-5.2830446
+1.5096,0.25,0.6,-4.57915
+1.5096,0.25,0.5,-3.846010
+1.5096,0.25,0.4,-3.145625
+1.5096,0.25,0.3,-2.446133
+1.5096,0.25,0.2,-1.718330
+1.5096,0.25,0.1,-0.883410
+1.5096,0.75,0.8,-7.528813
+1.5096,0.75,0.7,-6.89008
+1.5096,0.75,0.6,-6.0058145
+1.5096,0.75,0.55,-5.52616
+1.6438,0.5,0.8,-6.544010
+1.6438,0.5,0.7,-5.95657
+1.6438,0.5,0.6,-5.18836
+1.6438,0.5,0.5,-4.35859
+1.6438,0.5,0.4,-3.549019
+1.6438,0.5,0.3,-2.768327
+1.6438,0.5,0.2,-1.936931
+1.6438,0.5,0.1,-1.004111
diff --git a/validation/lj-mixture/run-validation.py b/validation/lj-mixture/run-validation.py
new file mode 100644
index 0000000..992e09e
--- /dev/null
+++ b/validation/lj-mixture/run-validation.py
@@ -0,0 +1,224 @@
+# /// script
+# requires-python = ">=3.13"
+# dependencies = [
+#     "matplotlib>=3.10.8",
+#     "numpy>=2.4.1",
+#     "pandas>=2.3.3",
+#     "seaborn>=0.13.2",
+# ]
+# ///
+
+# Comes from Monte Carlo Simulations of Binary Lennard–Jones Mixtures: A Test of the van der Waals One-Fluid Model, https://doi.org/10.1023/A:1022614200488
+
+import itertools
+import json
+import math
+import os
+import subprocess
+
+import matplotlib.pyplot as plt
+import numpy as np
+import pandas as pd
+import seaborn as sns
+
+# Different from publication, because we start from rectangular lattice instead of fcc
+# It shouldn't affect results, as we're not looking at crystallisation
+N = 1000
+path_to_exe = "particlesmc"
+
+epsilon_1 = 1
+epsilon_12 = 1.1523
+epsilon_2 = 1.3702
+sigma_1 = 1
+sigma_12 = 1.0339
+sigma_2 = 1.0640
+
+
+def create_config(n1: int, n2: int, box_length: float, filepath: str) -> None:
+    # We do a square lattice
+    # Species are put randomly amongst the sites
+
+    with open(filepath, "w") as f:
+        f.write(f"{N:d}\n")
+        f.write(
+            f'Lattice="{box_length:.4f} 0.0 0.0 0.0 {box_length:.4f} 0.0 0.0 0.0 {box_length:.4f}" Properties=:species:S:1:pos:R:3\n'
+        )
+
+        species_indices = [1] * n1 + [2] * n2
+        np.random.shuffle(species_indices)
+
+        number_of_particle_in_each_direction = round(N ** (1 / 3))
+        if number_of_particle_in_each_direction**3 != N:
+            raise RuntimeError("N is not x^3")
+        dxdydz = box_length / number_of_particle_in_each_direction
+        counter = 0
+        for i, j, k in itertools.product(
+            range(number_of_particle_in_each_direction),
+            range(number_of_particle_in_each_direction),
+            range(number_of_particle_in_each_direction),
+        ):
+            x = (i + 0.5) * dxdydz - box_length / 2
+            y = (j + 0.5) * dxdydz - box_length / 2
+            z = (k + 0.5) * dxdydz - box_length / 2
+            f.write(f"{species_indices[counter]} {x} {y} {z}\n")
+            counter += 1
+
+
+def create_params(parameters: dict, path_to_params: str) -> None:
+    with open(path_to_params, "w") as f:
+        f.write(f"""
+        [system]
+        config = "{parameters["config"]}"
+        temperature = {parameters["temperature"]}
+        density = {parameters["density"]}
+        list_type = "LinkedList"
+
+        [model]
+        [model."1-1"]
+        name = "LennardJones"
+        epsilon = {parameters["epsilon_1"]:.2f}
+        sigma = {parameters["sigma_1"]:.2f}
+        rcut = {parameters["rcut"]:.2f}
+        shift_potential = false
+
+        [model."1-2"]
+        name = "LennardJones"
+        epsilon = {parameters["epsilon_12"]}
+        sigma = {parameters["sigma_12"]}
+        rcut = {parameters["rcut"]:.2f}
+        shift_potential = false
+
+        [model."2-2"]
+        name = "LennardJones"
+        epsilon = {parameters["epsilon_2"]}
+        sigma = {parameters["sigma_2"]}
+        rcut = {parameters["rcut"]:.2f}
+        shift_potential = false
+
+        [simulation]
+        type = "Metropolis"
+        steps = {int(parameters["steps"]):d}
+        seed = 42
+        parallel = false
+        output_path = "./"
+
+        [[simulation.move]]
+        action = "Displacement"
+        probability = 1.0
+        policy = "SimpleGaussian"
+        parameters = {{sigma = 0.05}}
+
+        [[simulation.output]]
+        algorithm = "StoreCallbacks"
+        callbacks = ["energy", "acceptance"]
+        scheduler_params = {{linear_interval = 100}}
+
+        [[simulation.output]]
+        algorithm = "StoreLastFrames"
+        scheduler_params = {{linear_interval = 1000}}
+        fmt = "EXYZ"
+
+        """)
+
+
+def run_simulations(output_path: str) -> None:
+    df = pd.read_csv("reference-data.csv")
+
+    path_to_config = "config.exyz"
+    path_to_params = "params.toml"
+    data = []
+    for i, row in df.iterrows():
+        workdir = f"./tmp/{i}"
+        os.makedirs(workdir, exist_ok=True)
+
+        print(row["t"], row["x"], row["density"])
+
+        # density = N / box_length**3
+        box_length = (N / row["density"]) ** (1 / 3)
+
+        n2 = round(N * row["x"])
+        n1 = N - n2
+        create_config(n1, n2, box_length, f"{workdir}/{path_to_config}")
+
+        # Generate input, and run
+        rc = 4 * sigma_1
+        parameters = {
+            "config": path_to_config,
+            "temperature": row["t"],
+            "epsilon_1": epsilon_1,
+            "epsilon_12": epsilon_12,
+            "epsilon_2": epsilon_2,
+            "sigma_1": sigma_1,
+            "sigma_12": sigma_12,
+            "sigma_2": sigma_2,
+            "steps": 1e3,  # 1e4 equilibration in paper
+            "rcut": rc,
+            "density": N / box_length**3,
+        }
+        create_params(parameters, f"{workdir}/{path_to_params}")
+
+        subprocess.run(
+            [path_to_exe, path_to_params], cwd=workdir, stdout=subprocess.DEVNULL
+        )
+
+        # Post-process the energies
+        energies = pd.read_csv(f"{workdir}/energy.dat", sep="\\s+", names=["i", "e"])[
+            "e"
+        ]
+        # Remove the first half as equilibration, just to be sure
+        energies = energies[int(len(energies) / 2) :]
+
+        acceptance_rate = np.array(
+            pd.read_csv(f"{workdir}/acceptance.dat", sep="\\s+", names=["i", "a"])["a"]
+        )[-1]
+
+        # Compute long-range corrections from the cutoff
+        # Formula from Gromacs https://manual.gromacs.org/current/reference-manual/functions/long-range-vdw.html
+        lr_correction = 0
+
+        c6_11 = 4 * epsilon_1 * sigma_1**6
+        rho_11 = n1 / box_length**3
+        lr_correction += -2 / 3 * math.pi * n1 * rho_11 * c6_11 / rc**3
+
+        c6_22 = 4 * epsilon_2 * sigma_2**6
+        rho_22 = n2 / box_length**3
+        lr_correction += -2 / 3 * math.pi * n2 * rho_22 * c6_22 / rc**3
+
+        c6_12 = 4 * epsilon_12 * sigma_12**6
+        lr_correction += -2 / 3 * math.pi * (n1 * rho_22 + n2 * rho_11) * c6_12 / rc**3
+
+        data.append(
+            {
+                "t": row["t"],
+                "x": row["x"],
+                "density": row["density"],
+                "energy": np.mean(energies),
+                "energy_err": np.std(energies) / np.sqrt(len(energies)),
+                "lr_correction": lr_correction,
+                "acceptance_rate": json.loads(acceptance_rate)[0],
+            }
+        )
+
+        df.merge(pd.DataFrame(data)).to_csv(output_path)
+
+
+def plot_results(path_to_energies: str) -> None:
+    df = pd.read_csv(path_to_energies)
+
+    # u (from the paper) is actually total energy / epsilon_11 N
+    # epsilon_11 = 1, so it's the energy per particule
+    # Add long range corrections to the MC results, which is also in energy / particle
+    df["energy_lr"] = df["energy"] + df["lr_correction"] / N
+
+    sns.scatterplot(data=df, x="u", y="energy_lr", hue="density", style="x")
+    plt.axline((-5, -5), slope=1)
+    plt.xlabel("Published energy")
+    plt.ylabel("Re-computed from ParticleMC")
+    plt.savefig("correlation-plot.jpeg")
+    plt.show()
+
+
+if __name__ == "__main__":
+    path_to_energies = "calculated-energies.csv"
+    run_simulations(path_to_energies)
+    plot_results(path_to_energies)

From f4c063076a0b9c68120df5736605e50014cee16e Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Fran=C3=A7ois=20Villemot?= <francois.villemot@espci.fr>
Date: Thu, 15 Jan 2026 12:01:51 +0100
Subject: [PATCH 2/3] docs: README

---
 validation/lj-mixture/README.md | 30 ++++++++++++++++++++++++++++++
 1 file changed, 30 insertions(+)
 create mode 100644 validation/lj-mixture/README.md

diff --git a/validation/lj-mixture/README.md b/validation/lj-mixture/README.md
new file mode 100644
index 0000000..02d7024
--- /dev/null
+++ b/validation/lj-mixture/README.md
@@ -0,0 +1,30 @@
+# Binary mixture of LJ particles
+
+This directory contains scripts to run a simple simulation of a binary mixture of Lennard-Jones particles, in 3D.
+
+The average energy is recorded, for different conditions (temperature, density, ratio between the two species), and compared to published results.
+
+Reproducing the published results validates the implementation of the MC scheme, the LJ potential, and interaction counting.
+
+## How to run
+
+The `particlemc` exectuable should be in the PATH env variable, otherwise a specific path can
+be specified at the top of the run-validation.py script.
+
+The `run-validation.py` script is meant to be run with [uv](https://docs.astral.sh/uv/), with the following command:
+```bash
+uv run --script run-validation.py
+```
+which automatically install all dependencies.
+
+The script will run simulations at all desired parameters, save the results to csv, and make a correlation plot with published results.
+
+## Results
+
+ParticleMC reproduces published results very well:
+![correlation-plot](correlation-plot.jpeg)
+
+For some setups, the energies from ParticleMC are a bit higher than published ones. This actually happens at low density, because the simulation times are too short and we don't reach proper convergence (timeseries of energies not shown).
+
+## Reference data
+Comes from Monte Carlo Simulations of Binary Lennard–Jones Mixtures: A Test of the van der Waals One-Fluid Model, https://doi.org/10.1023/A:1022614200488

From 2c000dda69c2e5c40f1187d3bb7c76d90e36994d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Fran=C3=A7ois=20Villemot?= <francois.villemot@espci.fr>
Date: Thu, 15 Jan 2026 12:03:59 +0100
Subject: [PATCH 3/3] docs: format

---
 validation/lj-mixture/README.md | 4 ++++
 1 file changed, 4 insertions(+)

diff --git a/validation/lj-mixture/README.md b/validation/lj-mixture/README.md
index 02d7024..82dea55 100644
--- a/validation/lj-mixture/README.md
+++ b/validation/lj-mixture/README.md
@@ -22,8 +22,12 @@ The script will run simulations at all desired parameters, save the results to c
 ## Results
 
 ParticleMC reproduces published results very well:
+
+
 ![correlation-plot](correlation-plot.jpeg)
 
+On this plot, `x` is the ratio between the two types of particles, and the energy is per particle.
+
 For some setups, the energies from ParticleMC are a bit higher than published ones. This actually happens at low density, because the simulation times are too short and we don't reach proper convergence (timeseries of energies not shown).
 
 ## Reference data