From 086e7f7902c27585e80addb3338e6f85832f9013 Mon Sep 17 00:00:00 2001 From: Massimiliano Bonomi <massimiliano.bonomi@gmail.com> Date: Tue, 16 May 2017 11:33:08 +0200 Subject: [PATCH] finishing exercise 3 --- user-doc/figs/trieste-4-block-phi.png | Bin 0 -> 6042 bytes user-doc/tutorials/a-trieste-4.txt | 52 ++++++++++++++++++++++++-- 2 files changed, 49 insertions(+), 3 deletions(-) create mode 100644 user-doc/figs/trieste-4-block-phi.png diff --git a/user-doc/figs/trieste-4-block-phi.png b/user-doc/figs/trieste-4-block-phi.png new file mode 100644 index 0000000000000000000000000000000000000000..0125e2bf4d0b8806016574610591dd2d7d354953 GIT binary patch literal 6042 zcmb7Ic|26#`=2{T24QSjB12J(P?C_XgbCS~?2@JIJ6Yyt$&wIRvJ6s`ecv;KXtAZm znvAl`7}*Q+yZU^7-_P%l?{9m}%xli`zV|%mKJVw;=bZP6)w`y_#K6k{gTa`vnra3x z7#xJb01}c0Lf|aoLy#e(r)_kFN~PM^*i24N5(orkWn}>Yf!y5OwzjrBDmsL!(nqXp zr?%(S_fe>nz6Ly%N~BUjDv?Od^Q7k4Qu9=(pb7}WKwBzkOFU?*aDW`p7iR1tnWU#O zkvH)O_K_S;rapt^f#{8mFOByYUpBZ<X5$yq?)(D1saIeW5TxSts2nM<{y{1X`cR|$ zmte3y7!2g$@r1$d!(dd2`=&zog<w=1>?UNvsKc-dGMQT6PM!!SQK=MTO4-P45;=)- z_AC_6%*;%0Z!eKZtf{Ga{P-~rhjVmv)X~wAmX_w_<%PmmA-$(kVb8?GaMVvwG%Dc1 z7*e(%5v0OE6)F*}dTZ*<8yHaq1|yS+(0CAj@CVDw%Y(W@Lv#=1ltL}Q*Wj8VOl5Y& z5i;mKHLv@?V9c)$en2XQ`6vv=xsO#-H44mL&Y`a_5JRbq$NIJIH8OvpvdJ@_VMe`B z6BXLRuTgtoj|Rt2F@GX3t31SwF~Wi(z&r+cxEK;yj7C2X|My0SP~I(9tCwF*{C5q< z3l*L6(E>`S_hyFm^&xF#nX1HhnN!LzdI97JM<Vre*LwemK5;7IiJV(Inq?ZERzJB{ zpK?3Z%|ZdkfL4>RmH83cSIZ*N<BqDV=A-GL3&D^gREe+gBkLN`bK;ygJA>exFKv3x zpx<m;5_cl7#A1tEayJBn`;T}kB7+#r?jSg1srvjKyRI3!*E-h~#&^#3o*)H!I|lHI zJP^QU_;!Zigzm6WR$tholDiX6XgAxnHF3yZzrZkp8{@X~J*nIxLBFP&df!8Rds?HT zOhrYH7Et9xAdFEcBh~+t;kzCDUJ&hk_YYmiaMDO5!@t-&fB2G0iZTMB^OrY=pDchX z9U^idFw@6C8X0JD?d;4GM^?(_O-!EGOdKYf{>74xKB6@Vg9x!|zfkaOe`UsSy}#+( zRgLS2*0?q8rstze2;aFl>3$ADWFTica(gSlyaty}Ki|(e%sKwjgA2v4BvK82QyKmP zeXhfp^1-dUy5P-i4`n|S8qJ7gIjy#`9#NZNwgZoQ7XG+7<b~E4Vbc1~NNucVFWu|g z(GJwsLt!7rttCmj?U#tRvt%xyubYKs`^5+db(Iiz;0XeOWBf`985sBJN0wd<bLW7T zmmK?MMGO}-aFVtIv_*E7P2alPR{mVH0I0jg+c)H{DirUZ)>5D;7xSR?l@-94T22_9 zoQZ3ZZC{K$H~7+HZlHZp%XG&bljjnYv)s-0J_l1S-v5p>z(a!d{b+_3!4R=t*HXG6 zHqeJ$3_M>{2$CH#CM(+L**y&z0QL?H={2*&Z$1yLBILZe8DjU$(D)7lm@5v;VnM79 zU#Gjk;i-tWX^=`&Co-a**Gz<jIo+Jn<?!5?deThTOZQ<0WM4-h1A$z|`D#_Zo+<qp zklECPHJsM6RyuvOHV#bxX7IY;qkE0G&ArIcbe}h<6L}GzlPi=)(ku>*sQJ{zfkof0 zPwt=1X<2`QFQ5$=fk7T!?>A07j2jo+aC|(8oc|R&epL|}2ruDmv+684d(6tjB6xN6 zw{2G0&h8CkpIiEKY?O_BL$>sD66xprmemK!>*-Jwgo(=+;jXX(`8b98pXFY<DJ6~# zoJSmnlkvhNM}Bd0>*#<43M<8;-o<z0Yn@vsBD>+c2}iC(iR3P0lb+qqOk0M9rfG1w zx@9{nh0Wi>>Esnknq}x#$imQfP*!Wl{AmEeK|)%#vABxr{4(P!XDwqQI9WZuD1_bJ z{Goi<ZQEnTMDww`bm>mkWSLSJr!QlpOEC2)ilWW@N-k;o!4bMTQ(WtV)>}&xEMunL zeCA|C$R1|3JkKf+v5Vgj3tj0mhefOXc2un6_v8O_-@vMoygB_uJ~&Z-8gY1P`CWkX zRElHC*qs*Fj#iJ`<8OlXFTFz;$&0maFInKm>rqmX4VcNPsD{_q(-!>UZMR1HZtpN( z9F;H~sjwbq5M35)vKs0N{^QD~Y-?F$+?73m#|e?^UEbGve}T&6(bI?<9eVuzGYC3; z(|c~wIyP2ou(+Y@6tQythjrVGQ0Y?Cn~f;YO7j7jTPr5}@(S_t5XMA=RGn9WtpdLu z!&dJl(tg>UuR(k{cMps{ae}mY)FU;?hJC-Klh!v;HGBEm3W|bCr&l<0>5+FskY4e) z!dkhzGo|2cmQ+6c-tHa5BQ_!1xr-@mAEs4bNB}?t*Lh81T^};AcL&Z;zp!b{r9<rh zElISw@QSeE&I=?q8K|@iwFcu@1?aio&ptuJkNuR*+89#yN<of9fMIDRxG-E?xJ_e= zd2L{17Uu^x_#v+xFem)d<Xw<KJTmy@isi5-YbE0m&yU=2cdt#Karo#42t-!UBi5x_ zxK7yyXe27B{Y+#96t`)2#~)oHp3WsRp*-)p%*BCrXi`IKl}Cfgpu5H;xO@ud=DDtV zXL)$rFbi;b+qvx3Rg;kopN;)IkPJ<I2{zXEIi-~=^eOC~b{VD$mW(1yfZ|(*-|I|o zm&2Wg{mFv1;x~jbKT>zGz?AuonYbl|?pnlQUN_?BS8SEar*;8=GySBYkT8ZXfu0{O zCW3jT<R>J8DP7h)`<7v0pH0x9x??pf(!YINS_1i8t@LJ;Vl1-p*zYHJiK=AI&4YZ@ zGCDt;>O=BvK2x)XpCxN1xFReCb_8w3rp0YrNhDzmu(1E(I2~XGJ&e$_82SVpTJ0LR zP+K2FFo4`oEf%Qm|71fYfBy>qX8jYY;@%d(oV^l4JJW(OCUgN>9^Wa_PJdY&0HGE? z?W;BcIG~zp!T9D`!+KK==o?1p>O<UI`R($>j94>%4+E_xGfMF4tv&2t^?ybL?!CNh z+LRAE)zP16i~@ae1g%<2ujOAK^XBu>q!l)(Dvldjwd54+`T_&Dn1vu9TV~=Mt_ghy z<r}7~(`xAzDusw-8uac`wB}0x5k9z56q{P?!X4$Rpn`}Dgex@@t?>8QX)opr>!A(G z8|Ql5ipwRvyO@Rs_x8@R`gA?WkoM^kB%O|VgmM2=?R?&0fQ@qBQ6=51>DyxA55mJM zTb+G$0fOe`jjz#-XFBG?P!x&X`935vkgJ2XbnJ=5*u);fRDCM~3=q^KR!>F@`R#2} z$*w>B1}CLrqeBVTGC*U)tc9`67f0DK&AG~mR(d;lIeEHd(QGiZ)&C`OduewX9eJGa z8zo}fg1+G^upJ=Kfv9B^QV-}D*)rMvD!~nSM5z&r<+sGq@NgBP{|JzK2N?sR#Ph1$ zj7`d9SPNR(+Z0&EZSC(VI|gffP92;;lX8$p-n;L2$B-lUR?3hoz-vJB*;kq;lV!*w zm2TRjsk5JgS(Q*p`zwx-iQ)I{I9fx?mp{xiQ&t?2{WNKXGUF$EDjtwwfQNb`=pe&( z=GT#Hz*m;`AeB9rc4VM(herDqISAz?w-b&}%>lYfyNt1IhkzB$O7QoBPDNtdq2FYv z?O>R3H)fFQ{oceQ3j1Jl1PakAIu(($86H8@MO?&*s<5v1A<q;@st~i-l#pju3nbgy zf`u-yS~)9mYHS4<b_BF%MxZtc4KuR%xXXlbHj%u+SfvPDgRnNSduWQiCdj8tMx8q} z5}@SJ=ivJ|bPN#_%+A>d>HFmyXN_o;4!0*3m$~_TVU_pk3YEsn!N2lMp`VV))y`Oz zH~JrC0%|mbY<ax7=La~4SE3yTPgbsUM>1}6Ta1LZ_3(Qb_!$qqd7_IZJ$PtIf2Ljm zpV=u%QZmkV=HfHp2xA`n#CM+5?(0$BSUsw{eAM?A?H)Pqlxeor?iRC6`2`Is*C(0} z)QO|FSOPD~J(M=h9w&GLi9PKeLXF3pLw{J&6wXTZTY4qzt$3{;L$u0%pcjK6&zZPj zC1Y{z;di3M@#C49ov-kjG=PWh*ezZZW%FrHA#h;GKD{xyj5SkfCh#F~Vq|RXJTil( z9eza7&&T~q5(*fIr+X_nLEBK&Ke&{5Z7ZGL0O--InLchS5<?sB%X<ir7Th(b#H``N z8J|EP(n;Mp&e?GNhRu-^!N-!wk}N0S>U~`CU=aExKgPVxcDIM7{W4a6jrO~ZNR0d2 zm2WeQ>O^gn?%>h2i=@Z-Xq}E=V&hKusw8`du9`+4mo{An<5U4DV<%jWE`{kybn=sf z4z8MOSCRV;movb%<@a4p!<VuLn^*)-M#9ezSsvpDd@eOJCU^?*LXD8afqSveh+g_t zX~8r2xX!kaB>EAv;Kfr2C}$$1(Hi!qf*NXHj4e=vuZX8;LtfRkC!`|o7RGwY9)>%r zZTBKsXnMwmhLs;h;^YLIc=9WC-lJaSlzHw4a~upMg%oui5#|2gYjT0A7@q*j2|gCA zt5TJIDqq9dk4r50ZgVYL`Es#>#gx3uRLi_Cc-?D%*^KNOwNx6l0&!+ydD?dC-Dge_ zv{U|?5QQ7sKK<?_N+UuDdY;ES5aOAe3x+}GRMTWF^0kP2s>1Nu2=sM3<MHOe!*}m* zh$?Zlxkd>pq62O?Xl&P!qJxLdFR#9E4yeC`KyYzgd~`E=lz|mxX}M~Yjj^zc21!D7 zug_LVN<;M%tF5w0Xkt}Toq1d1-3UFXd*meR1JD!rume;)NWhJ<Dk<lkE+5E907nJ2 zjeNYieG$D0Pu1JE(*>}?#j3<zO;LTOlYsPzXKFIr>4KsM*kP=QK9g+6H1=%lZIqO5 zGS(;p3`*91!i(r-Z=aAdJz3n|#)qAE=SK8Cj3nPvCGuhC+_`sql6?;h@dE?Vo9v6H zZzqUjojBsb>&d;RG88@WC!`tFi5Ibs5EKMK+|-|HLbgG_?9gyW-A1r5=tr1jut=n$ zXCV?#FB-or5iEiLq_3=N$D~fMF2CiOgHrQY=#I$rjtzHrj&wf{Kh^un7!3nH&Sx&( zM?-^X$(kMGj==<pyn?Q>Lz%rIFMVlST9BOWRX%@_Ir&WZ8u;u?vc-*7(H8$fRw?`P z7PfE@Mcjb*3j{(7fNdetyDJowb^78>Uq~29Kh5<up47uL7Xz{9oey@WNcXey12r@6 z#;Y#Medilom9+TKD1ZWLAbFH}0x^w}viWOG6XpRbm=G`iR*Ry@H-fg^f9apAVsPuD zbjkndDXfZ0H(Ks*kr|5&bmP|K`b%jIawT0Rmj11`9Dv8^rvm=+elXg(LxRKVcO(B% zV-hIOR3Z1UkLB#I`v!ULPWSD=GiQ2R*RsB?c8rwWMP?|@F*aUh-3l>ej=gvPjsQQf zHtP$Rj`GaPex`J|7F=freC8m4K6zHQr55W0&66xmc-oNf-d)lx8jA5`GH5kuv=@_c z3^imaLfX_B9JxVYtUm&*%_LiTkyYQkk(iKIwB5Ay<aK-uJ&C<e!&N?6C@Zb8{zA>7 z0*X`TCuQJ{15xSU4z12=c_rEt@1%P`bZYE5A`rbtp`dGNrwu%3|0aq7B8n#Y(Elbp z4Knh6j&UcQ`%9lEgW&Kzu=XEX{>9jTBVUEht1o@mHi-e7*j*ClmGtdebWLJ#P3)rA zr?ypVvV;WypOYr|q|~$usH!Sd8Tf6q*lM)dB6@So;Qmo{m-|V65L1vb6k6Spb3pUM zWKNSI-_XTOxAka{P7Agre51OfN~6?&^y~uU4BD{0RWoHdl?am+vaaqZF5yi8>@L^d zjB>LN*cU`U{~Aw;21oy%4~uYXURW~Szc?6VuOuAJ<uR!F2eL!K|2~+A^Dn>fhVTEQ z(*KK?7n+GH`SHd-VV)&3+Je<3?knCqw^-}z#V_}vJ<y#vb#8vLZmFf0mFx_|L$!dP zs((%w)(rd<<vXnCGT!_-M7YLLS*y-F@phnGt0P2jMhytKemkK+VQ172o550?or)>< zxmiLkylIyt<HR_BG_rA%2&02w_phQfnip2oYL`$b61e8wwLEZ%hR++?62MmeO(s*d zXwUToZzVqe(cEu0Q(yJqq4C0!Yur8Y?}4cbtykqNY09qn6^}8SrgoHw<HM0GHIhO5 zSJ&LJrga#or<RC@(F3c`a;28C#x|M;GgV?G?pEh+@>p2Jq`5X<xolKKc-+9YpNXFm zjiJbQZ({P;qC}4g$LA#kXvvpYiINum_h{nXtQ$V9`z44aHmI(r+k}_+=EzJ0?wIQp zt=wXGo4i5&P{cLmcW)(<P?%M$*?BJN>D61m6Pli1UfOPWS~C^}FL%$9%d1}YJRMZd zZnmS3HGdXds;8e*xICQw^eUlFez1avF{b$N;he&}{*Tc6e*2W&CDU4*R&$BfDbiy9 zXPS53t&i=XzD@2=jo*FKS238Y>i*Ry^1@*FJ75)<^V%fMgy*<S&hu*o%rceTXp}3| z;x$2d;xij-<&qMc=SA;iqaZ%Lp7TD`6%(<%vL21QCVF#KlV>|w%YF>UZ`od5<Ve)N z?@ODJZ_R0WQUK6;!833*Tz(*e>sCfs<#E<Z<@dLr?d=3IW={_vJ}&!K`N6(a&Eo5d zDT0U$Ky$z1h{AL4(N7=5?DNjHUi*FV4OjFF>(Fj5erbMZ`n1v|{9JU$Qanq4y~(Ca zT=$XWD{HRUfI789f)lkP$@Y|6#}2Kcr>i77mp%=oxDDSE&4tHx#0P7%^K!T3E~eP{ zT*wT3_zrzDbY)-gJNah0k693&PEBvdZ&9SNAC-@adZC2;$nj0AV&b}I&RYP(VBc3i zugX7}=4+xgtF~S{;*-<)PUj_)4_}4*yD2@Jz~%;4e9S(HFUZES_2xhkx?^emwU4o6 zI9hwv?3SD1z2Q|itW%w@L|0w4me*_wqrQMu+i#MAr67yFl3DeRC~|Z|Ym}2P#a$^u z>(^>=cJtW<@bdpweEwMl{&Sk|hEM6&ciGIPYvejG!^7^Hsf^qOZ9s4`0<^Kc@1^u$ z;f42jMtPP%%e=C9p{UKTAlY9<g2!@o2C6fzM3t%%B`fPsKaeSK`It@9$>F&d;gmS6 zOTD*~^y_2{{^bqjXId^rbEKb#0D^5sM73M@PQA9Q&X7Z)<oGDzw%`p}juLzg-XV8_ z4!XhfdG$=C6aO6F^VzqeZuWCdm05*lvkU{Vr;3_TLXtkW(IhUjk4zsFTX}P%PP@rq z9NZMgH$wCRQI*+v9|8I>_g{s7Zh5%?j1NyU?V0W|2car+@ps!+=%LthM0{8IXs)H` z>#W-xj(Ak{F$;4OidhwRzH6lKR3mP1#o_m`y_0Orh2GxgN8TBQQ7KW#ZQA-WH~X*I xbmB!`F)zw;`e}}Q&$p?{C=_=2AKBhvS+YZ(3;upM6FMG-VXs_ME5BqL{vRsh=kx#o literal 0 HcmV?d00001 diff --git a/user-doc/tutorials/a-trieste-4.txt b/user-doc/tutorials/a-trieste-4.txt index 7a47600f9..21ac0feac 100644 --- a/user-doc/tutorials/a-trieste-4.txt +++ b/user-doc/tutorials/a-trieste-4.txt @@ -393,12 +393,58 @@ dihedral \f$ \phi \f$ along with the corresponding weight: 1.130800 0.0622998 \endverbatim +At this point we can apply the block-analysis technique we have learned in the +\ref trieste-2 tutorial to calculate for different block sizes the average free-energy +and the error. For your convenience, you can use the `do_block_fes_norm.py` python +script to read the `phi.weight` file and produce the desired output. +We use a bash loop to use block sizes ranging from 1 to 1000: + +\verbatim +for i in `seq 1 10 1000`; do python3 do_block_fes.py phi.weight 1 -3.141593 3.018393 51 2.494339 $i; done +\endverbatim + +For each value of block length `N`, you will obtain a separate `fes.N.dat` file, containing the value +of the \f$ \phi \f$ variable on a grid, the average free-energy, and the associated error (in Kjoule/mol): + +\verbatim + -3.141593 23.184653 0.080659 + -3.018393 17.264462 0.055181 + -2.895194 13.360259 0.047751 + -2.771994 10.772696 0.043548 + -2.648794 9.403544 0.042022 +\endverbatim + +Finally, we can calculate the average error along the free-energy profile as a function of the block length: + +\verbatim +for i in `seq 1 10 1000`; do ave=`awk '{tot+=$3}END{print tot/NR}' fes.$i.dat`; echo $i $ave; done > err.blocks +\endverbatim + +and visualize it using `gnuplot`: + +\verbatim +gnuplot> p "err.blocks" u 1:2 w lp +\endverbatim + +As expected, the error increases with the block length until it reaches a plateau in correspondence of a dimension +of the block that exceeds the correlation between data points (Fig. \ref trieste-4-block-phi). + +\anchor trieste-4-block-phi +\image html trieste-4-block-phi.png "Block analysis of a metadynamics simulation using phi as CV" + +To finish this exercize, calculate the error associated to the free energy as a function of +the collective variable \f$ \psi \f$ computed in Exercize 2 and compare with the error along \f$ \phi \f$. +What can we learn from this analysis about the convergence of the two metadynamics simulations +and the quality of the collective variables chosen? + \section trieste-4-conclusions Conclusions In summary, in this tutorial you should have learned how to use PLUMED to: -- Manipulate atomic coordinates. -- Compute collective variables. - +- Setup and run a metadynamics calculation. +- Compute free energies from the metadynamics bias potential using the \ref sum_hills utility. +- Calculate the error in the reconstructed free energy using block analysis. +- Discriminate between good and bad collective variables. +- Evaluate the convergence of a metadynamics simulation. */ -- GitLab