From db2040af09e73f3fc871158afe11ee242ec6cfb1 Mon Sep 17 00:00:00 2001 From: Adam Sieradzan Date: Wed, 30 Jan 2019 16:28:14 +0100 Subject: [PATCH] changes in wham and cluser + cutoff corr --- PARAM/contact_ext.3.parm | 576 ++++++++++++++++++++++++++ source/cluster/CMakeLists.txt | 3 + source/cluster/cluster.F90 | 51 ++- source/cluster/io_clust.F90 | 8 +- source/cluster/main_clust.F | 41 +- source/unres/data/control_data.F90 | 2 +- source/unres/energy.F90 | 46 ++- source/wham/CMakeLists.txt | 3 + source/wham/conform_compar.F90 | 8 +- source/wham/enecalc.F90 | 22 + source/wham/io_wham.F90 | 779 ++++++++++++++++++++++++++++-------- source/wham/wham_calc.F90 | 16 +- 12 files changed, 1335 insertions(+), 220 deletions(-) diff --git a/PARAM/contact_ext.3.parm b/PARAM/contact_ext.3.parm index 8e59ad1..784188a 100644 --- a/PARAM/contact_ext.3.parm +++ b/PARAM/contact_ext.3.parm @@ -574,3 +574,579 @@ 0.35531356 -0.07100358 3.61044593 2.4 ala dbz 0.32491593 0.00199416 3.48648915 1.4 ala aib 0.32491593 0.00199416 3.48648915 1.4 ala abu + 0.35860954 0.16982996 3.29863666 3.2 cys cys + 0.56282835 -0.05299321 3.60729199 2.1 cys met + 0.33964241 0.05207525 4.17039421 2.4 cys phe + 0.37589252 -0.13754074 3.40782436 3.2 cys ile + 0.19324876 -0.12584486 3.97277369 2.4 cys leu + 0.40072817 0.10036272 4.49587991 1.6 cys val + 0.42263839 0.06477987 3.63916207 3.6 cys trp + 0.54237237 0.04744217 3.48308370 3.1 cys tyr + 0.51547365 -0.00799778 3.40011763 2.2 cys ala + 0.51607945 -0.12446741 3.70421858 1.6 cys gly + 0.43412348 0.15581630 4.56079263 1.3 cys thr + 0.38622584 -0.07788669 3.77359853 2.0 cys ser + 0.50860385 -0.01971684 3.49518888 3.0 cys gln + 0.39057873 -0.13582425 3.09416287 3.6 cys asn + 0.42554013 -0.10224673 3.99873956 2.0 cys glu + 0.51197275 -0.00799778 3.42940555 3.0 cys asp + 0.51368036 -0.01971684 3.45279388 2.8 cys his + 0.53687925 -0.05416896 3.53908504 2.8 cys arg + 0.52071055 -0.00432467 3.85035238 2.3 cys lys + 0.56790413 -0.10114534 3.62096650 1.9 cys pro + 0.56282835 -0.05299321 3.60729199 2.1 cys sem + 0.33964241 0.05207525 4.17039421 2.4 cys phe + 0.51547365 -0.00799778 3.40011763 2.2 cys aib + 0.51547365 -0.00799778 3.40011763 2.2 cys abu + 0.68048894 -0.06899195 3.60729199 2.1 met cys + 0.47168906 -0.01389221 3.62448254 3.4 met met + 0.38809483 0.05839704 3.69650927 3.8 met phe + 0.30403358 -0.17995008 4.29852296 2.5 met ile + 0.33973020 -0.34834021 4.94950826 1.2 met leu + 0.30071237 0.01285368 4.18347745 2.6 met val + 0.52333332 0.17569346 3.51061576 4.2 met trp + 0.61170525 0.20950877 4.06080648 3.0 met tyr + 0.41902197 -0.08043824 3.17380260 3.3 met ala + 0.52139540 0.00438353 3.43153041 2.6 met gly + 0.49784550 0.12355497 4.19683352 2.2 met thr + 0.63180941 -0.16149377 4.09465313 1.6 met ser + 0.44886809 -0.15991724 3.53814108 3.3 met gln + 0.41073744 0.00740397 4.06717512 2.5 met asn + 0.54442625 -0.07267571 3.41282057 3.5 met glu + 0.63996137 -0.06915184 3.94820588 2.2 met asp + 0.62501187 0.11051556 3.67829283 2.8 met his + 0.53066821 0.05034669 3.92004157 3.0 met arg + 0.33560623 -0.08136851 4.49980555 2.0 met lys + 0.49366704 -0.10767906 3.76955189 2.7 met pro + 0.47168906 -0.01389221 3.62448254 3.4 met sem + 0.38809483 0.05839704 3.69650927 3.8 met dbz + 0.41902197 -0.08043824 3.17380260 3.3 met aib + 0.41902197 -0.08043824 3.17380260 3.3 met abu + 0.36262117 -0.19622749 4.17039421 2.4 phe cys + 0.26128467 0.01513555 3.69650927 3.8 phe met + 0.20490480 -0.09214195 4.22979013 3.0 phe phe + 0.28320764 -0.11272094 4.28523432 2.8 phe ile + 0.49776022 0.15233796 3.52640364 3.7 phe leu + 0.38397299 0.13538890 4.83243966 2.0 phe val + 0.37298546 -0.07130682 3.38863637 4.4 phe trp + 0.50532840 -0.04912285 3.56784740 3.7 phe tyr + 0.47086276 -0.20131132 3.61044593 2.4 phe ala + 0.47136787 -0.02932642 3.89454305 2.0 phe gly + 0.43636452 0.05389401 4.20612459 2.6 phe thr + 0.39664510 0.04040127 3.97909758 2.6 phe ser + 0.37270815 -0.11268458 3.88510036 3.0 phe gln + 0.55039146 -0.03205935 4.10110530 2.6 phe asn + 0.33746013 -0.24962819 3.93827274 2.7 phe glu + 0.43948791 -0.14279533 4.16332864 2.3 phe asp + 0.57384992 -0.03571288 3.97249455 2.6 phe his + 0.49181945 0.18810757 3.92791573 3.4 phe arg + 0.34273176 -0.21712089 3.95551322 3.0 phe lys + 0.44953604 -0.00243978 3.77053096 2.9 phe pro + 0.26128467 0.01513555 3.69650927 3.8 phe sme + 0.20490480 -0.09214195 4.22979013 3.0 phe dbz + 0.47086276 -0.20131132 3.61044593 2.4 phe aib + 0.47086276 -0.20131132 3.61044593 2.4 phe abu + 0.40831404 -0.11390349 3.40782436 3.2 ile cys + 0.36166967 0.16294633 4.29852296 2.5 ile met + 0.17153929 0.01989158 4.28523432 2.8 ile phe + 0.47535754 0.30766747 3.94168510 3.2 ile ile + 0.29127023 -0.18034538 3.95620992 2.9 ile leu + 0.40891029 0.37221078 4.77230778 2.0 ile val + 0.40490198 0.25444624 3.55185603 4.3 ile trp + 0.40189859 0.12737180 3.99860644 3.1 ile tyr + 0.37442735 -0.07448267 3.66929044 2.7 ile ala + 0.38837277 0.21898895 3.87773858 2.1 ile gly + 0.48188659 0.47591913 4.59886714 2.2 ile thr + 0.26566057 -0.00988237 4.43023881 1.8 ile ser + 0.16126661 0.12105441 3.21763584 4.2 ile gln + 0.37751545 -0.07665185 3.55147963 3.1 ile asn + 0.26828456 -0.05165526 3.82344372 3.1 ile glu + 0.11152774 -0.21162156 4.20071886 2.1 ile asp + 0.27999892 -0.04581434 3.87438904 3.0 ile his + 0.28625871 0.18884379 4.27655600 2.9 ile arg + 0.51810565 0.20449041 4.52497861 1.8 ile lys + 0.27899493 0.14745809 4.93089469 1.5 ile pro + 0.36166967 0.16294633 4.29852296 2.5 ile sme + 0.17153929 0.01989158 4.28523432 2.8 ile dbz + 0.37442735 -0.07448267 3.66929044 2.7 ile aib + 0.37442735 -0.07448267 3.66929044 2.7 ile abu + 0.27669006 -0.08857650 3.97277369 2.4 leu cys + 0.46145819 0.11704558 4.94950826 1.2 leu met + 0.39499219 -0.06849736 3.52640364 3.7 leu phe + 0.37951422 0.12929207 3.95620992 2.9 leu ile + 0.37548431 0.25972538 4.68988538 2.2 leu leu + 0.31811044 -0.12698715 4.07045515 2.6 leu val + 0.56177403 0.27868668 3.78215175 3.7 leu trp + 0.38385722 0.09687124 3.56586168 3.7 leu tyr + 0.34823804 -0.09467021 4.02083643 1.8 leu ala + 0.47111001 0.19515789 3.78987702 2.1 leu gly + 0.48994191 0.30052632 4.61644076 2.0 leu thr + 0.33626266 -0.13893079 3.78032608 2.5 leu ser + 0.47005465 -0.09199330 3.89621932 2.5 leu gln + 0.41288855 0.02311019 3.71732693 2.8 leu asn + 0.32910008 -0.16223127 4.40927943 1.8 leu glu + 0.29408977 -0.05874025 4.12366306 2.3 leu asp + 0.43376746 0.01891589 3.48566489 3.5 leu his + 0.44700953 0.11595371 3.53775960 3.6 leu arg + 0.38429980 0.03881280 4.36103401 2.2 leu lys + 0.43262307 0.03816639 4.08852337 2.4 leu pro + 0.46145819 0.11704558 4.94950826 1.2 leu sme + 0.39499219 -0.06849736 3.52640364 3.7 leu dbz + 0.34823804 -0.09467021 4.02083643 1.8 leu aib + 0.34823804 -0.09467021 4.02083643 1.8 leu abu + 0.26123252 0.11055566 4.49587991 1.6 val cys + -0.11325692 -0.14606457 4.18347745 2.6 val met + 0.09409692 0.11177909 4.83243966 2.0 val phe + -0.00666800 -0.07497911 4.77230778 2.0 val ile + 0.28991366 0.17674772 4.07045515 2.6 val leu + -0.10515666 -0.11003497 4.12314563 2.8 val val + 0.30183179 0.14806296 3.74692371 3.7 val trp + 0.29338214 0.08282836 4.24592089 2.6 val tyr + 0.12064673 -0.00982537 4.10663351 1.7 val ala + 0.19385817 -0.03823541 3.80654563 2.0 val gly + -0.01013325 -0.07812848 3.79210798 3.1 val thr + 0.29858743 0.07767388 3.58093493 2.8 val ser + 0.35822386 -0.09616652 3.47357789 3.2 val gln + 0.34627205 -0.01311854 3.51658941 3.1 val asn + 0.06254100 -0.04565289 4.06702052 2.5 val glu + 0.18901220 -0.16654191 4.12654645 2.1 val asp + 0.24391223 0.13836662 4.17589586 2.4 val his + 0.06426380 0.09368312 4.00988606 3.0 val arg + -0.13595797 -0.17161034 4.64018900 2.0 val lys + 0.17932755 0.09383245 4.56454004 2.0 val pro + -0.11325692 -0.14606457 4.18347745 2.6 val sme + 0.09409692 0.11177909 4.83243966 2.0 val dbz + 0.12064673 -0.00982537 4.10663351 1.7 val aib + 0.12064673 -0.00982537 4.10663351 1.7 val abu + 0.31196785 0.02064366 3.63916207 3.6 trp cys + 0.46366596 -0.03700266 3.51061576 4.2 trp met + 0.58088610 -0.11066174 3.38863637 4.4 trp phe + 0.24583477 -0.33494927 3.55185603 4.3 trp ile + 0.32671974 -0.22702523 3.78215175 3.7 trp leu + 0.45612392 0.10634243 3.74692371 3.7 trp val + 0.36417512 -0.20528168 2.98506467 5.8 trp trp + 0.50155918 -0.12768864 3.15117399 5.1 trp tyr + 0.32812767 -0.00610991 3.47226781 3.7 trp ala + 0.40467084 -0.29644951 3.34936384 3.4 trp gly + 0.49007720 -0.23396191 4.06972159 2.9 trp thr + 0.46206287 -0.16415403 3.19708779 4.0 trp ser + 0.42889946 -0.18850961 3.34772285 4.3 trp gln + 0.51128091 -0.12547768 3.16744534 4.5 trp asn + 0.39266635 -0.37772802 3.09551785 4.5 trp glu + 0.48243867 0.01467772 3.56336594 3.9 trp asp + 0.42702047 -0.16605199 3.41608597 4.4 trp his + 0.52626838 -0.12547768 3.01681706 5.0 trp arg + 0.37799932 -0.12350334 3.76243382 3.8 trp lys + 0.49906417 -0.08578986 3.94460593 2.7 trp pro + 0.46366596 -0.03700266 3.51061576 4.2 trp sme + 0.58088610 -0.11066174 3.38863637 4.4 trp dbz + 0.32812767 -0.00610991 3.47226781 3.7 trp aib + 0.32812767 -0.00610991 3.47226781 3.7 trp abu + 0.61872403 -0.14590233 3.48308370 3.1 tyr cys + 0.51373838 0.10165708 4.06080648 3.0 tyr met + 0.48847918 -0.25213997 3.56784740 3.7 tyr phe + 0.61169095 0.09973185 3.99860644 3.1 tyr ile + 0.59378113 -0.01121051 3.56586168 3.7 tyr leu + 0.46928231 -0.00248010 4.24592089 2.6 tyr val + 0.43672853 -0.13447255 3.15117399 5.1 tyr trp + 0.48884849 -0.22303211 2.84305042 5.2 tyr tyr + 0.67005094 -0.07584168 3.22022624 3.4 tyr ala + 0.65174425 -0.24311887 3.86011481 2.0 tyr gly + 0.50879141 -0.25935026 4.15837027 2.5 tyr thr + 0.54564025 -0.22727470 3.66401120 2.7 tyr ser + 0.66853458 -0.14863634 3.50128149 3.5 tyr gln + 0.65954068 -0.14062835 3.33886112 3.6 tyr asn + 0.50689915 -0.32717288 3.42836329 3.5 tyr glu + 0.51002994 -0.14640234 3.25938544 4.0 tyr asp + 0.61200946 -0.24411836 4.02269765 2.4 tyr his + 0.46534272 -0.20073886 3.36685862 4.3 tyr arg + 0.42351858 -0.18191847 4.14326778 2.9 tyr lys + 0.47769600 -0.16559256 3.95285757 2.9 tyr pro + 0.51373838 0.10165708 4.06080648 3.0 tyr sme + 0.48847918 -0.25213997 3.56784740 3.7 tyr dbz + 0.67005094 -0.07584168 3.22022624 3.4 tyr aib + 0.67005094 -0.07584168 3.22022624 3.4 tyr abu + 0.31170564 0.13087300 3.40011763 2.2 ala cys + 0.35817150 -0.05118436 3.17380260 3.3 ala met + 0.35531356 -0.07100358 3.61044593 2.4 ala phe + 0.16898670 0.06047415 3.66929044 2.7 ala ile + 0.29707520 -0.00251129 4.02083643 1.8 ala leu + 0.16009076 -0.19966586 4.10663351 1.7 ala val + 0.33543361 0.02729744 3.47226781 3.7 ala trp + 0.40771637 0.04944451 3.22022624 3.4 ala tyr + 0.32491593 0.00199416 3.48648915 1.4 ala ala + 0.32250163 0.00199416 3.50576712 1.1 ala gly + 0.44809712 -0.03625897 3.28319535 2.9 ala thr + 0.32570984 0.00199416 3.48018928 2.2 ala ser + 0.25708601 -0.02518611 3.98829105 2.0 ala gln + 0.27223871 0.02170282 3.77047757 2.1 ala asn + 0.37126794 -0.02710099 3.50057256 2.6 ala glu + 0.15955468 -0.05289865 3.56513816 2.4 ala asp + 0.41120950 0.04944451 3.19268981 2.9 ala his + 0.27446733 0.00786488 4.09477647 1.9 ala arg + 0.22981425 -0.09322596 4.42291744 0.9 ala lys + 0.37507388 0.05588781 3.38612759 2.5 ala pro + 0.35817150 -0.05118436 3.17380260 3.3 ala sme + 0.35531356 -0.07100358 3.61044593 2.4 ala dbz + 0.32491593 0.00199416 3.48648915 1.4 ala aib + 0.32491593 0.00199416 3.48648915 1.4 ala abu + -0.00336656 0.02622444 3.70421858 1.6 gly cys + -0.09302725 0.26319430 3.43153041 2.6 gly met + -0.09760057 0.17790940 3.89454305 2.0 gly phe + -0.22751502 0.17812915 3.87773858 2.1 gly ile + -0.25069250 0.09450769 3.78987702 2.1 gly leu + -0.18213085 0.07188032 3.80654563 2.0 gly val + 0.02178246 0.08557004 3.34936384 3.4 gly trp + -0.11493048 0.22893358 3.86011481 2.0 gly tyr + -0.11323243 0.10823838 3.50576712 1.1 gly ala + -0.11186659 0.10823838 3.52493966 0.8 gly gly + -0.15040905 0.28755373 3.59967226 2.3 gly thr + -0.11368436 0.10823838 3.49950195 1.9 gly ser + -0.03597803 0.10164763 4.26660952 1.4 gly gln + -0.11204539 0.10823838 3.52240926 2.4 gly asn + -0.23143037 -0.03794759 3.90940249 1.9 gly glu + -0.14158106 -0.07444570 3.87388976 1.8 gly asp + -0.11031945 0.23711722 3.23796583 3.0 gly his + -0.17946054 0.00419438 3.25922115 3.2 gly arg + -0.29030826 0.10161339 4.01153241 1.5 gly lys + -0.01947913 0.15568872 3.41316198 2.1 gly pro + -0.09302725 0.26319430 3.43153041 2.6 gly sme + -0.09760057 0.17790940 3.89454305 2.0 gly dbz + -0.11323243 0.10823838 3.50576712 1.1 gly aib + -0.11323243 0.10823838 3.50576712 1.1 gly abu + 0.25327931 0.01671477 4.56079263 1.3 thr cys + 0.23416404 0.05555785 4.19683352 2.2 thr met + -0.05802816 -0.16538168 4.20612459 2.6 thr phe + -0.00442492 -0.07151185 4.59886714 2.2 thr ile + 0.07792594 0.14574366 4.61644076 2.0 thr leu + 0.17777770 0.04200606 3.79210798 3.1 thr val + 0.28036614 0.17671473 4.06972159 2.9 thr trp + 0.21670162 0.02128130 4.15837027 2.5 thr tyr + 0.28703884 0.05657710 3.28319535 2.9 thr ala + 0.20185014 -0.02758099 3.59967226 2.3 thr gly + -0.02168311 0.07585571 4.09544630 2.8 thr thr + 0.12417560 -0.13819171 3.43029089 3.1 thr ser + 0.19988109 -0.12547841 3.29954085 3.6 thr gln + 0.09516162 -0.02303978 3.87932811 2.9 thr asn + -0.00509545 -0.16916448 3.91402729 2.8 thr glu + -0.23145136 -0.34804909 3.89588321 2.3 thr asp + 0.01936770 -0.15186953 4.11716100 2.7 thr his + 0.16868560 0.15199778 4.12871646 2.8 thr arg + 0.37511648 0.23392103 5.19110983 0.8 thr lys + 0.14543985 0.10143310 3.81482913 2.9 thr pro + 0.23416404 0.05555785 4.19683352 2.2 thr sme + -0.05802816 -0.16538168 4.20612459 2.6 thr dbz + 0.28703884 0.05657710 3.28319535 2.9 thr aib + 0.28703884 0.05657710 3.28319535 2.9 thr abu + 0.21517295 -0.24554423 3.77359853 2.0 ser cys + 0.25883244 -0.12670009 4.09465313 1.6 ser met + 0.14102918 -0.01538099 3.97909758 2.6 ser phe + 0.25921833 0.08059688 4.43023881 1.8 ser ile + 0.20136778 -0.16976156 3.78032608 2.5 ser leu + 0.28531354 -0.09273778 3.58093493 2.8 ser val + 0.39448436 -0.03696454 3.19708779 4.0 ser trp + 0.43603784 0.09648460 3.66401120 2.7 ser tyr + 0.30692270 -0.04308303 3.48018928 2.2 ser ala + 0.30457335 -0.04308303 3.49950195 1.9 ser gly + 0.32020038 0.00223073 3.43029089 3.1 ser thr + 0.30769547 -0.04308303 3.47387798 2.4 ser ser + 0.11488573 -0.13472978 3.12392468 3.9 ser gln + 0.21105760 -0.12158253 3.72665420 2.4 ser asn + 0.04115496 0.04092516 4.61451625 1.5 ser glu + 0.10088016 -0.28980003 3.72754795 2.4 ser asp + 0.17721632 -0.10937848 3.18728475 3.7 ser his + 0.13355810 -0.32587585 3.87794402 2.5 ser arg + 0.22605173 -0.17882693 4.37211026 1.5 ser lys + 0.13546239 -0.16030706 4.22414116 2.0 ser pro + 0.25883244 -0.12670009 4.09465313 1.6 ser sme + 0.14102918 -0.01538099 3.97909758 2.6 ser dbz + 0.30692270 -0.04308303 3.48018928 2.2 ser aib + 0.30692270 -0.04308303 3.48018928 2.2 ser abu + 0.46907450 0.23814067 3.49518888 3.0 gln cys + 0.64462489 0.20629501 3.53814108 3.3 gln met + 0.62450221 0.26002714 3.88510036 3.0 gln phe + 0.64883105 0.31118595 3.21763584 4.2 gln ile + 0.52131440 0.03272724 3.89621932 2.5 gln leu + 0.65106047 0.16613237 3.47357789 3.2 gln val + 0.57768605 0.14932176 3.34772285 4.3 gln trp + 0.55788468 0.30256567 3.50128149 3.5 gln tyr + 0.47680255 0.04342046 3.98829105 2.0 gln ala + 0.47868436 -0.03204437 4.26660952 1.4 gln gly + 0.45754127 0.06306618 3.29954085 3.6 gln thr + 0.59818612 0.19256561 3.12392468 3.9 gln ser + 0.55149012 -0.02069686 3.92298950 2.4 gln gln + 0.55970840 0.15394656 3.33545255 3.7 gln asn + 0.52123754 0.23540666 3.53114474 3.5 gln glu + 0.56477134 0.12111678 4.15342990 2.3 gln asp + 0.60162954 0.19256561 3.09986815 4.1 gln his + 0.41096673 -0.00430890 3.59510761 3.4 gln arg + 0.54772483 0.12647211 4.40523407 1.8 gln lys + 0.65056464 0.45016870 4.27763332 2.2 gln pro + 0.64462489 0.20629501 3.53814108 3.3 gln sme + 0.62450221 0.26002714 3.88510036 3.0 gln dbz + 0.47680255 0.04342046 3.98829105 2.0 gln aib + 0.47680255 0.04342046 3.98829105 2.0 gln abu + 0.57572212 0.14373398 3.09416287 3.6 asn cys + 0.44333739 0.28309684 4.06717512 2.5 asn met + 0.44442019 0.21261136 4.10110530 2.6 asn phe + 0.49499497 0.11891284 3.55147963 3.1 asn ile + 0.53576677 0.15090652 3.71732693 2.8 asn leu + 0.52751327 0.11911315 3.51658941 3.1 asn val + 0.52684562 0.21607903 3.16744534 4.5 asn trp + 0.57007352 0.24209638 3.33886112 3.6 asn tyr + 0.53360274 0.10469691 3.77047757 2.1 asn ala + 0.53096441 0.17749158 3.52240926 2.4 asn gly + 0.40760213 0.22182945 3.87932811 2.9 asn thr + 0.47722200 0.11590346 3.72665420 2.4 asn ser + 0.43939620 0.14592084 3.33545255 3.7 asn gln + 0.46547953 0.08929379 3.07503648 4.2 asn asn + 0.45066424 0.03916198 3.88381092 2.6 asn glu + 0.56591130 0.16761648 3.73937544 2.6 asn asp + 0.53142890 0.30777945 3.38670300 3.4 asn his + 0.30442390 0.01364112 4.48141529 1.9 asn arg + 0.27593110 0.18004212 4.18933053 2.5 asn lys + 0.50353752 0.32164896 3.73467727 2.9 asn pro + 0.44333739 0.28309684 4.06717512 2.5 asn sme + 0.44442019 0.21261136 4.10110530 2.6 asn dbz + 0.53360274 0.10469691 3.77047757 2.1 asn aib + 0.53360274 0.10469691 3.77047757 2.1 asn abu + 0.45167162 0.28718023 3.99873956 2.0 glu cys + 0.58329544 0.31669128 3.41282057 3.5 glu met + 0.65649724 0.21234517 3.93827274 2.7 glu phe + 0.53672118 0.18971451 3.82344372 3.1 glu ile + 0.61529042 0.23787159 4.40927943 1.8 glu leu + 0.57711221 0.16371932 4.06702052 2.5 glu val + 0.57842612 0.14725703 3.09551785 4.5 glu trp + 0.68045805 0.20875844 3.42836329 3.5 glu tyr + 0.52828434 0.09169166 3.50057256 2.6 glu ala + 0.54349648 0.27262886 3.90940249 1.9 glu gly + 0.65080974 0.43402569 3.91402729 2.8 glu thr + 0.47242896 0.07323090 4.61451625 1.5 glu ser + 0.45962325 0.34422613 3.53114474 3.5 glu gln + 0.45784628 0.08542837 3.88381092 2.6 glu asn + 0.45136983 0.32035022 3.77778958 2.9 glu glu + 0.40995790 -0.03116378 4.26798075 2.0 glu asp + 0.52859120 0.16354546 3.78696151 3.1 glu his + 0.55740234 0.22045560 3.81136201 3.1 glu arg + 0.55181260 0.35744091 4.66520100 1.6 glu lys + 0.66236078 0.36052823 3.47970840 3.3 glu pro + 0.58329544 0.31669128 3.41282057 3.5 glu sme + 0.65649724 0.21234517 3.93827274 2.7 glu dbz + 0.52828434 0.09169166 3.50057256 2.6 glu aib + 0.52828434 0.09169166 3.50057256 2.6 glu abu + 0.32481843 0.30241943 3.42940555 3.0 asp cys + 0.48146154 0.12184204 3.94820588 2.2 asp met + 0.49531096 0.13459205 4.16332864 2.3 asp phe + 0.48629858 0.10892335 4.20071886 2.1 asp ile + 0.42458969 0.04067298 4.12366306 2.3 asp leu + 0.46012374 0.06388031 4.12654645 2.1 asp val + 0.44400377 0.21493328 3.56336594 3.9 asp trp + 0.49469528 0.08857244 3.25938544 4.0 asp tyr + 0.43791837 0.09152666 3.56513816 2.4 asp ala + 0.40782681 0.10505865 3.87388976 1.8 asp gly + 0.53260595 0.06977694 3.89588321 2.3 asp thr + 0.41899901 0.06697791 3.72754795 2.4 asp ser + 0.31470323 0.16090621 4.15342990 2.3 asp gln + 0.41087617 0.03408678 3.73937544 2.6 asp asn + 0.36889970 -0.15930847 4.26798075 2.0 asp glu + 0.49914534 0.28237257 3.82849221 2.6 asp asp + 0.35831041 0.17803592 3.53238450 3.3 asp his + 0.46844169 0.00518452 2.71889129 5.3 asp arg + 0.33974035 0.00746795 4.11707156 2.3 asp lys + 0.58065436 0.12893457 3.46466683 2.9 asp pro + 0.48146154 0.12184204 3.94820588 2.2 asp sme + 0.49531096 0.13459205 4.16332864 2.3 asp dbz + 0.43791837 0.09152666 3.56513816 2.4 asp aib + 0.43791837 0.09152666 3.56513816 2.4 asp abu + 0.66712991 -0.00059312 3.45279388 2.8 his cys + 0.62023666 0.07810136 3.67829283 2.8 his met + 0.54813450 0.02387399 3.97249455 2.6 his phe + 0.50448776 0.10129703 3.87438904 3.0 his ile + 0.43502367 -0.13271534 3.48566489 3.5 his leu + 0.50083039 0.09390639 4.17589586 2.4 his val + 0.34876524 0.17810403 3.41608597 4.4 his trp + 0.52214146 -0.08932064 4.02269765 2.4 his tyr + 0.62795477 -0.02594442 3.19268981 2.9 his ala + 0.55372190 0.05824969 3.23796583 3.0 his gly + 0.41891546 0.11138432 4.11716100 2.7 his thr + 0.54068082 -0.04376077 3.18728475 3.7 his ser + 0.52837481 -0.09497389 3.09986815 4.1 his gln + 0.59926199 -0.09807584 3.38670300 3.4 his asn + 0.30784556 -0.28464757 3.78696151 3.1 his glu + 0.57435087 0.02111146 3.53238450 3.3 his asp + 0.65809056 -0.09152492 2.98047276 4.0 his his + 0.54935375 -0.06444900 3.84130645 3.2 his arg + 0.55747594 0.18348617 4.00467347 2.6 his lys + 0.61284214 0.16197710 3.89481022 2.5 his pro + 0.62023666 0.07810136 3.67829283 2.8 his sme + 0.54813450 0.02387399 3.97249455 2.6 his dbz + 0.62795477 -0.02594442 3.19268981 2.9 his aib + 0.62795477 -0.02594442 3.19268981 2.9 his abu + 0.75225350 0.20735450 3.53908504 2.8 arg cys + 0.61621105 0.32748378 3.92004157 3.0 arg met + 0.66251194 0.14629924 3.92791573 3.4 arg phe + 0.57629387 0.00597023 4.27655600 2.9 arg ile + 0.74108830 0.24297317 3.53775960 3.6 arg leu + 0.62859689 0.04738249 4.00988606 3.0 arg val + 0.76107791 0.21620637 3.01681706 5.0 arg trp + 0.64334219 -0.01420347 3.36685862 4.3 arg tyr + 0.61765505 -0.23266856 4.09477647 1.9 arg ala + 0.63577289 0.03826170 3.25922115 3.2 arg gly + 0.52962288 0.04658043 4.12871646 2.8 arg thr + 0.57322475 -0.01630823 3.87794402 2.5 arg ser + 0.63768385 -0.04885203 3.59510761 3.4 arg gln + 0.60938603 0.03300984 4.48141529 1.9 arg asn + 0.63351724 -0.10980006 3.81136201 3.1 arg glu + 0.63619125 -0.04461328 2.71889129 5.3 arg asp + 0.61095825 0.10189057 3.84130645 3.2 arg his + 0.71730775 0.08728438 2.81437659 5.4 arg arg + 0.60658407 -0.03333816 4.23716346 2.5 arg lys + 0.68568576 0.07786298 3.77385723 3.0 arg pro + 0.61621105 0.32748378 3.92004157 3.0 arg sme + 0.66251194 0.14629924 3.92791573 3.4 arg dbz + 0.61765505 -0.23266856 4.09477647 1.9 arg aib + 0.61765505 -0.23266856 4.09477647 1.9 arg abu + 0.76470157 0.76084085 3.85035238 2.3 lys cys + 0.67255968 0.44875322 4.49980555 2.0 lys met + 0.76960024 0.60831176 3.95551322 3.0 lys phe + 0.61378093 0.26242515 4.52497861 1.8 lys ile + 0.71563444 0.68072580 4.36103401 2.2 lys leu + 0.67409525 0.58814366 4.64018900 2.0 lys val + 0.88671210 1.00380282 3.76243382 3.8 lys trp + 0.67158296 0.42635325 4.14326778 2.9 lys tyr + 0.67254432 0.35248523 4.42291744 0.9 lys ala + 0.69462202 0.31735268 4.01153241 1.5 lys gly + 0.58760205 0.22591520 5.19110983 0.8 lys thr + 0.62437605 0.35173136 4.37211026 1.5 lys ser + 0.70291910 0.54587675 4.40523407 1.8 lys gln + 0.71487329 0.62683782 4.18933053 2.5 lys asn + 0.63139916 0.33371478 4.66520100 1.6 lys glu + 0.72656677 0.46183912 4.11707156 2.3 lys asp + 0.83364969 0.89847717 4.00467347 2.6 lys his + 0.77862679 0.71215250 4.23716346 2.5 lys arg + 0.56582949 -0.04350298 4.78186650 1.3 lys lys + 0.80330046 0.84741487 4.53882609 2.0 lys pro + 0.67255968 0.44875322 4.49980555 2.0 lys sme + 0.76960024 0.60831176 3.95551322 3.0 lys dbz + 0.67254432 0.35248523 4.42291744 0.9 lys aib + 0.67254432 0.35248523 4.42291744 0.9 lys abu + 0.49304685 -0.20581200 3.62096650 1.9 pro cys + 0.12683702 -0.30939905 3.76955189 2.7 pro met + 0.11511445 -0.18829660 3.77053096 2.9 pro phe + 0.05125111 -0.18929142 4.93089469 1.5 pro ile + 0.11747631 -0.02693818 4.08852337 2.4 pro leu + 0.26952062 0.30671521 4.56454004 2.0 pro val + 0.38994700 -0.04061465 3.94460593 2.7 pro trp + 0.29377104 0.11826092 3.95285757 2.9 pro tyr + 0.25804259 -0.13915496 3.38612759 2.5 pro ala + 0.38094914 -0.14696919 3.41316198 2.1 pro gly + 0.30752420 0.17160109 3.81482913 2.9 pro thr + 0.07822398 0.05027273 4.22414116 2.0 pro ser + 0.14107113 -0.04129291 4.27763332 2.2 pro gln + 0.21471187 -0.00268809 3.73467727 2.9 pro asn + 0.22170637 -0.11224658 3.47970840 3.3 pro glu + 0.08872349 -0.40993001 3.46466683 2.9 pro asp + 0.31135187 0.03163507 3.89481022 2.5 pro his + 0.38918420 0.21697911 3.77385723 3.0 pro arg + -0.11536833 -0.12070488 4.53882609 2.0 pro lys + 0.24887782 -0.07176013 3.66640730 2.7 pro pro + 0.12683702 -0.30939905 3.76955189 2.7 pro sme + 0.11511445 -0.18829660 3.77053096 2.9 pro dbz + 0.25804259 -0.13915496 3.38612759 2.5 pro aib + 0.25804259 -0.13915496 3.38612759 2.5 pro abu + 0.68048894 -0.06899195 3.60729199 2.1 sme cys + 0.47168906 -0.01389221 3.62448254 3.4 sme met + 0.38809483 0.05839704 3.69650927 3.8 sme phe + 0.30403358 -0.17995008 4.29852296 2.5 sme ile + 0.33973020 -0.34834021 4.94950826 1.2 sme leu + 0.30071237 0.01285368 4.18347745 2.6 sme val + 0.52333332 0.17569346 3.51061576 4.2 sme trp + 0.61170525 0.20950877 4.06080648 3.0 sme tyr + 0.41902197 -0.08043824 3.17380260 3.3 sme ala + 0.52139540 0.00438353 3.43153041 2.6 sme gly + 0.49784550 0.12355497 4.19683352 2.2 sme thr + 0.63180941 -0.16149377 4.09465313 1.6 sme ser + 0.44886809 -0.15991724 3.53814108 3.3 sme gln + 0.41073744 0.00740397 4.06717512 2.5 sme asn + 0.54442625 -0.07267571 3.41282057 3.5 sme glu + 0.63996137 -0.06915184 3.94820588 2.2 sme asp + 0.62501187 0.11051556 3.67829283 2.8 sme his + 0.53066821 0.05034669 3.92004157 3.0 sme arg + 0.33560623 -0.08136851 4.49980555 2.0 sme lys + 0.49366704 -0.10767906 3.76955189 2.7 sme pro + 0.47168906 -0.01389221 3.62448254 3.4 sme sme + 0.38809483 0.05839704 3.69650927 3.8 sme dbz + 0.41902197 -0.08043824 3.17380260 3.3 sme aib + 0.41902197 -0.08043824 3.17380260 3.3 sme abu + 0.36262117 -0.19622749 4.17039421 2.4 dbz cys + 0.26128467 0.01513555 3.69650927 3.8 dbz met + 0.20490480 -0.09214195 4.22979013 3.0 dbz phe + 0.28320764 -0.11272094 4.28523432 2.8 dbz ile + 0.49776022 0.15233796 3.52640364 3.7 dbz leu + 0.38397299 0.13538890 4.83243966 2.0 dbz val + 0.37298546 -0.07130682 3.38863637 4.4 dbz trp + 0.50532840 -0.04912285 3.56784740 3.7 dbz tyr + 0.47086276 -0.20131132 3.61044593 2.4 dbz ala + 0.47136787 -0.02932642 3.89454305 2.0 dbz gly + 0.43636452 0.05389401 4.20612459 2.6 dbz thr + 0.39664510 0.04040127 3.97909758 2.6 dbz ser + 0.37270815 -0.11268458 3.88510036 3.0 dbz gln + 0.55039146 -0.03205935 4.10110530 2.6 dbz asn + 0.33746013 -0.24962819 3.93827274 2.7 dbz glu + 0.43948791 -0.14279533 4.16332864 2.3 dbz asp + 0.57384992 -0.03571288 3.97249455 2.6 dbz his + 0.49181945 0.18810757 3.92791573 3.4 dbz arg + 0.34273176 -0.21712089 3.95551322 3.0 dbz lys + 0.44953604 -0.00243978 3.77053096 2.9 dbz pro + 0.26128467 0.01513555 3.69650927 3.8 dbz sme + 0.20490480 -0.09214195 4.22979013 3.0 dbz dbz + 0.47086276 -0.20131132 3.61044593 2.4 dbz aib + 0.47086276 -0.20131132 3.61044593 2.4 dbz abu + 0.31170564 0.13087300 3.40011763 2.2 ala cys + 0.35817150 -0.05118436 3.17380260 3.3 ala met + 0.35531356 -0.07100358 3.61044593 2.4 ala phe + 0.16898670 0.06047415 3.66929044 2.7 ala ile + 0.29707520 -0.00251129 4.02083643 1.8 ala leu + 0.16009076 -0.19966586 4.10663351 1.7 ala val + 0.33543361 0.02729744 3.47226781 3.7 ala trp + 0.40771637 0.04944451 3.22022624 3.4 ala tyr + 0.32491593 0.00199416 3.48648915 1.4 ala ala + 0.32250163 0.00199416 3.50576712 1.1 ala gly + 0.44809712 -0.03625897 3.28319535 2.9 ala thr + 0.32570984 0.00199416 3.48018928 2.2 ala ser + 0.25708601 -0.02518611 3.98829105 2.0 ala gln + 0.27223871 0.02170282 3.77047757 2.1 ala asn + 0.37126794 -0.02710099 3.50057256 2.6 ala glu + 0.15955468 -0.05289865 3.56513816 2.4 ala asp + 0.41120950 0.04944451 3.19268981 2.9 ala his + 0.27446733 0.00786488 4.09477647 1.9 ala arg + 0.22981425 -0.09322596 4.42291744 0.9 ala lys + 0.37507388 0.05588781 3.38612759 2.5 ala pro + 0.35817150 -0.05118436 3.17380260 3.3 ala sme + 0.35531356 -0.07100358 3.61044593 2.4 ala dbz + 0.32491593 0.00199416 3.48648915 1.4 ala aib + 0.32491593 0.00199416 3.48648915 1.4 ala abu + 0.31170564 0.13087300 3.40011763 2.2 ala cys + 0.35817150 -0.05118436 3.17380260 3.3 ala met + 0.35531356 -0.07100358 3.61044593 2.4 ala phe + 0.16898670 0.06047415 3.66929044 2.7 ala ile + 0.29707520 -0.00251129 4.02083643 1.8 ala leu + 0.16009076 -0.19966586 4.10663351 1.7 ala val + 0.33543361 0.02729744 3.47226781 3.7 ala trp + 0.40771637 0.04944451 3.22022624 3.4 ala tyr + 0.32491593 0.00199416 3.48648915 1.4 ala ala + 0.32250163 0.00199416 3.50576712 1.1 ala gly + 0.44809712 -0.03625897 3.28319535 2.9 ala thr + 0.32570984 0.00199416 3.48018928 2.2 ala ser + 0.25708601 -0.02518611 3.98829105 2.0 ala gln + 0.27223871 0.02170282 3.77047757 2.1 ala asn + 0.37126794 -0.02710099 3.50057256 2.6 ala glu + 0.15955468 -0.05289865 3.56513816 2.4 ala asp + 0.41120950 0.04944451 3.19268981 2.9 ala his + 0.27446733 0.00786488 4.09477647 1.9 ala arg + 0.22981425 -0.09322596 4.42291744 0.9 ala lys + 0.37507388 0.05588781 3.38612759 2.5 ala pro + 0.35817150 -0.05118436 3.17380260 3.3 ala sme + 0.35531356 -0.07100358 3.61044593 2.4 ala dbz + 0.32491593 0.00199416 3.48648915 1.4 ala aib + 0.32491593 0.00199416 3.48648915 1.4 ala abu diff --git a/source/cluster/CMakeLists.txt b/source/cluster/CMakeLists.txt index cdb54a5..2c701e2 100644 --- a/source/cluster/CMakeLists.txt +++ b/source/cluster/CMakeLists.txt @@ -78,6 +78,9 @@ if(UNRES_MD_FF STREQUAL "GAB" ) elseif(UNRES_MD_FF STREQUAL "E0LL2Y") # set preprocesor flags set(CPPFLAGS "PROCOR -DSPLITELE -DSCCORPDB" ) +elseif(UNRES_MD_FF STREQUAL "NEWCORR") + # set preprocesor flags + set(CPPFLAGS "PROCOR -DUNRES -DISNAN -DSPLITELE -DLANG0 -DNEWCORR -DCORRCD" ) elseif(UNRES_MD_FF STREQUAL "4P") set(CPPFLAGS "SPLITELE -DLANG0 -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORPDB" ) endif(UNRES_MD_FF STREQUAL "GAB") diff --git a/source/cluster/cluster.F90 b/source/cluster/cluster.F90 index 3b75b28..70b685e 100644 --- a/source/cluster/cluster.F90 +++ b/source/cluster/cluster.F90 @@ -14,7 +14,7 @@ c,cref use energy_data, only: nnt,nct use control_data, only: symetr,outpdb,outmol2,titel,& - iopt,print_dist !,MaxProcs + iopt,print_dist,nclust !,MaxProcs use control, only: tcpu,initialize use wham_data, only: punch_dist @@ -54,7 +54,7 @@ INTEGER,dimension(:),allocatable :: HVALS !(maxconf-1) INTEGER,dimension(:),allocatable :: IORDER,HEIGHT !(maxconf-1) integer,dimension(:),allocatable :: nn !(maxconf) - integer :: ndis + integer :: ndis,is,ie real(kind=4),dimension(:),allocatable :: DISNN !(maxconf) LOGICAL,dimension(:),allocatable :: FLAG !(maxconf) integer :: i,j,k,l,m,n,len,lev,idum,ii,ind,jj,icut,ncon,& @@ -96,8 +96,9 @@ !elwrite(iout,*) "before parmread" call openunits !elwrite(iout,*) "before parmread" - call parmread call read_control + call parmread +! call read_control !elwrite(iout,*) "after read control" call molread ! if (refstr) call read_ref_structure(*30) @@ -114,7 +115,12 @@ ! write (iout,*) "after permut" ! call flush(iout) print *,'MAIN: nnt=',nnt,' nct=',nct - + if (nclust.gt.0) then + PRINTANG(1)=.TRUE. + PRINTPDB(1)=outpdb + printmol2(1)=outmol2 + ncut=0 + else DO I=1,NCUT PRINTANG(I)=.FALSE. PRINTPDB(I)=0 @@ -126,12 +132,21 @@ printmol2(i)=outmol2 ENDIF ENDDO + endif + if (ncut.gt.0) then write (iout,*) 'Number of cutoffs:',NCUT write (iout,*) 'Cutoff values:' DO ICUT=1,NCUT WRITE(IOUT,'(8HRCUTOFF(,I2,2H)=,F8.1,2i2)')ICUT,RCUTOFF(ICUT),& printpdb(icut),printmol2(icut) ENDDO + else if (nclust.gt.0) then + write (iout,'("Number of clusters requested",i5)') nclust + else + if (me.eq.Master) & + write (iout,*) "ERROR: Either nclust or ncut must be >0" + stop + endif DO I=1,NRES-3 MULT(I)=1 ENDDO @@ -327,21 +342,29 @@ allocate(iass(maxgr)) allocate(nconf(maxgr,maxingr)) allocate(totfree_gr(maxgr)) - +!c 3/3/16 AL: added explicit number of cluters + if (nclust.gt.0) then + is=nclust-1 + ie=nclust-1 + icut=1 + else + is=1 + ie=lev-1 + endif do i=1,maxgr licz(i)=0 enddo icut=1 - i=1 - NGR=i+1 + i=is + NGR=is+1 do j=1,n licz(iclass(j,i))=licz(iclass(j,i))+1 nconf(iclass(j,i),licz(iclass(j,i)))=j ! write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)), ! & nconf(iclass(j,i),licz(iclass(j,i))) enddo - do i=1,lev-1 - +! do i=1,lev-1 + do i=is,ie idum=lev-i DO L=1,LEV IF (HEIGHT(L).EQ.IDUM) GOTO 190 @@ -349,7 +372,8 @@ 190 IDUM=L write(IOUT,*) "i+1",i+1," idum",idum," critval",CRITVAL(IDUM),& " icut",icut," cutoff",rcutoff(icut) - IF (CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN + IF (nclust.gt.0.or.CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN + if (nclust.le.0) & WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut) write (iout,'(a,f8.2)') 'Maximum distance found:',& CRITVAL(IDUM) @@ -363,11 +387,10 @@ do l=1,maxgr licz(l)=0 enddo + ii=i-is+1 do j=1,n - enddo - do j=1,n - licz(iclass(j,i))=licz(iclass(j,i))+1 - nconf(iclass(j,i),licz(iclass(j,i)))=j + licz(iclass(j,ii))=licz(iclass(j,ii))+1 + nconf(iclass(j,ii),licz(iclass(j,ii)))=j !d write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)),& !d nconf(iclass(j,i),licz(iclass(j,i))) !d print *,j,iclass(j,i), diff --git a/source/cluster/io_clust.F90 b/source/cluster/io_clust.F90 index 1e24794..b4cf707 100644 --- a/source/cluster/io_clust.F90 +++ b/source/cluster/io_clust.F90 @@ -1257,7 +1257,7 @@ use control_data, only: titel,outpdb,outmol2,refstr,pdbref,& iscode,symetr,punch_dist,print_dist,nstart,nend,& caonly,iopt,efree,lprint_cart,lprint_int,rlamb_ele,& - r_cut_ele + r_cut_ele,nclust,tor_mode,scelemode ! implicit none ! include 'DIMENSIONS' ! include 'sizesclu.dat' @@ -1314,7 +1314,13 @@ call reada(controlcard,'BOXX',boxxsize,100.0d0) call reada(controlcard,'BOXY',boxysize,100.0d0) call reada(controlcard,'BOXZ',boxzsize,100.0d0) + call readi(controlcard,'NCLUST',nclust,5) ! ions=index(controlcard,"IONS").gt.0 + call readi(controlcard,"SCELEMODE",scelemode,0) + print *,"SCELE",scelemode + call readi(controlcard,'TORMODE',tor_mode,0) +!C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then + write(iout,*) "torsional and valence angle mode",tor_mode call reada(controlcard,"R_CUT_ELE",r_cut_ele,15.0d0) call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0) write(iout,*) "R_CUT_ELE=",r_cut_ele diff --git a/source/cluster/main_clust.F b/source/cluster/main_clust.F index 15e0bd0..be5e364 100644 --- a/source/cluster/main_clust.F +++ b/source/cluster/main_clust.F @@ -61,8 +61,9 @@ C call initialize call openunits - call parmread call read_control + call parmread +! call read_control call molread c if (refstr) call read_ref_structure(*30) do i=1,nres @@ -79,6 +80,12 @@ c write (iout,*) "after permut" c call flush(iout) print *,'MAIN: nnt=',nnt,' nct=',nct + if (nclust.gt.0) then + PRINTANG(1)=.TRUE. + PRINTPDB(1)=outpdb + printmol2(1)=outmol2 + ncut=0 + else DO I=1,NCUT PRINTANG(I)=.FALSE. PRINTPDB(I)=0 @@ -90,12 +97,21 @@ c call flush(iout) printmol2(i)=outmol2 ENDIF ENDDO + endif + if (ncut.gt.0) then write (iout,*) 'Number of cutoffs:',NCUT write (iout,*) 'Cutoff values:' DO ICUT=1,NCUT WRITE(IOUT,'(8HRCUTOFF(,I2,2H)=,F8.1,2i2)')ICUT,RCUTOFF(ICUT), & printpdb(icut),printmol2(icut) ENDDO + else if (nclust.gt.0) then + write (iout,'("Number of clusters requested",i5)') nclust + else + if (me.eq.Master) + & write (iout,*) "ERROR: Either nclust or ncut must be >0" + stop + endif DO I=1,NRES-3 MULT(I)=1 ENDDO @@ -241,19 +257,28 @@ C goto 192 endif CALL HCASS(N,IA,IB,CRIT,LEV,ICLASS,HVALS,IORDER,CRITVAL,HEIGHT) -c CALL HCDEN(LEV,IORDER,HEIGHT,CRITVAL) +!c CALL HCDEN(LEV,IORDER,HEIGHT,CRITVAL) +!c 3/3/16 AL: added explicit number of cluters + if (nclust.gt.0) then + is=nclust-1 + ie=nclust-1 + icut=1 + else + is=1 + ie=lev-1 + endif do i=1,maxgr licz(i)=0 enddo icut=1 - i=1 - NGR=i+1 + i=is + NGR=is+1 do j=1,n licz(iclass(j,i))=licz(iclass(j,i))+1 nconf(iclass(j,i),licz(iclass(j,i)))=j -c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)), -c & nconf(iclass(j,i),licz(iclass(j,i))) +!c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)), +!c & nconf(iclass(j,i),licz(iclass(j,i))) enddo do i=1,lev-1 @@ -264,8 +289,8 @@ c & nconf(iclass(j,i),licz(iclass(j,i))) 190 IDUM=L write(IOUT,*) "i+1",i+1," idum",idum," critval",CRITVAL(IDUM), & " icut",icut," cutoff",rcutoff(icut) - IF (CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN - WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut) + IF (nclust.gt.0.or.CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN + if (nclust.le.0) WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut) write (iout,'(a,f8.2)') 'Maximum distance found:', & CRITVAL(IDUM) CALL SRTCLUST(ICUT,ncon_work,iT) diff --git a/source/unres/data/control_data.F90 b/source/unres/data/control_data.F90 index e4fb4e9..4f2203f 100644 --- a/source/unres/data/control_data.F90 +++ b/source/unres/data/control_data.F90 @@ -37,7 +37,7 @@ gnorm_check,gradout,split_ene,with_theta_constr,protein,ions,nucleic,& raw_psipred #ifdef CLUSTER - integer :: iopt,nend,nstart,outpdb,outmol2 !cluster + integer :: iopt,nend,nstart,outpdb,outmol2,nclust !cluster logical :: punch_dist,print_dist,lside,lprint_cart,lprint_int,& caonly,efree,from_bx,from_cx,from_cart ! cluster #else diff --git a/source/unres/energy.F90 b/source/unres/energy.F90 index 1d24115..094c724 100644 --- a/source/unres/energy.F90 +++ b/source/unres/energy.F90 @@ -1,4 +1,4 @@ - module energy + module energy !----------------------------------------------------------------------------- use io_units use names @@ -793,6 +793,8 @@ call AFMforce(Eafmforce) else if (selfguide.gt.0) then call AFMvel(Eafmforce) + else + Eafmforce=0.0d0 endif endif if (tubemode.eq.1) then @@ -831,7 +833,7 @@ eespp=0.0d0 endif ! write(iout,*) ecorr_nucl,"ecorr_nucl",nres_molec(2) -! print *,"before ecatcat" +! print *,"before ecatcat",wcatcat if (nfgtasks.gt.1) then if (fg_rank.eq.0) then call ecatcat(ecationcation) @@ -916,8 +918,8 @@ ! Here are the energies showed per procesor if the are more processors ! per molecule then we sum it up in sum_energy subroutine ! print *," Processor",myrank," calls SUM_ENERGY" - energia(41)=ecation_prot - energia(42)=ecationcation + energia(42)=ecation_prot + energia(41)=ecationcation energia(46)=escbase energia(47)=epepbase energia(48)=escpho @@ -1042,8 +1044,8 @@ etors_d_nucl=energia(36) ecorr_nucl=energia(37) ecorr3_nucl=energia(38) - ecation_prot=energia(41) - ecationcation=energia(42) + ecation_prot=energia(42) + ecationcation=energia(41) escbase=energia(46) epepbase=energia(47) escpho=energia(48) @@ -1255,8 +1257,8 @@ etors_d_nucl=energia(36) ecorr_nucl=energia(37) ecorr3_nucl=energia(38) - ecation_prot=energia(41) - ecationcation=energia(42) + ecation_prot=energia(42) + ecationcation=energia(41) escbase=energia(46) epepbase=energia(47) escpho=energia(48) @@ -1331,11 +1333,11 @@ ecorr,wcorr,& ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3,& eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr,& - ethetacnstr,ebr*nss,Uconst,eliptran,wliptran,Eafmforc, & + ethetacnstr,ebr*nss,Uconst,eliptran,wliptran,Eafmforce, & etube,wtube, & estr_nucl,wbond_nucl, ebe_nucl,wang_nucl,& - evdwpp,wvdwpp_nucl,eespp,welpp,evdwpsb,wvdwpsb,eelpsb,welpsb& - evdwsb,wvdwsb,eelsb,welsb,esbloc,wsbloc,etors_nucl,wtor_nucl& + evdwpp,wvdwpp_nucl,eespp,welpp,evdwpsb,wvdwpsb,eelpsb,welpsb,& + evdwsb,wvdwsb,eelsb,welsb,esbloc,wsbloc,etors_nucl,wtor_nucl,& etors_d_nucl,wtor_d_nucl,ecorr_nucl,wcorr_nucl,& ecorr3_nucl,wcorr3_nucl,ecation_prot,wcatprot,ecationcation,wcatcat, & escbase,wscbase,epepbase,wpepbase,escpho,wscpho,epeppho,wpeppho,& @@ -2013,8 +2015,8 @@ ! write(iout,*)"c ", c(1,:), c(2,:), c(3,:) rrij=1.0D0/(xj*xj+yj*yj+zj*zj) rij=dsqrt(rrij) - sss_ele_cut=sscale_ele(1.0d0/(rij*sigma(itypi,itypj))) - sss_ele_grad=sscagrad_ele(1.0d0/(rij*sigma(itypi,itypj))) + sss_ele_cut=sscale_ele(1.0d0/(rij)) + sss_ele_grad=sscagrad_ele(1.0d0/(rij)) ! print *,sss_ele_cut,sss_ele_grad,& ! 1.0d0/(rij),r_cut_ele,rlamb_ele if (sss_ele_cut.le.0.0) cycle @@ -2076,7 +2078,7 @@ fac=rij*fac ! print *,'before fac',fac,rij,evdwij fac=fac+evdwij*sss_ele_grad/sss_ele_cut& - /sigma(itypi,itypj)*rij + *rij ! print *,'grad part scale',fac, & ! evdwij*sss_ele_grad/sss_ele_cut & ! /sigma(itypi,itypj)*rij @@ -11321,7 +11323,7 @@ +wcorr3_nucl*gradcorr3_nucl(j,i) +& wcatprot* gradpepcat(j,i)+ & wcatcat*gradcatcat(j,i)+ & - wscbase*gvdwc_scbase(j,i) & + wscbase*gvdwc_scbase(j,i)+ & wpepbase*gvdwc_pepbase(j,i)+& wscpho*gvdwc_scpho(j,i)+& wpeppho*gvdwc_peppho(j,i) @@ -11552,7 +11554,7 @@ +gradafm(j,i) & +wliptran*gliptranc(j,i) & +welec*gshieldc(j,i) & - +welec*gshieldc_loc(j,) & + +welec*gshieldc_loc(j,i) & +wcorr*gshieldc_ec(j,i) & +wcorr*gshieldc_loc_ec(j,i) & +wturn3*gshieldc_t3(j,i) & @@ -14020,8 +14022,8 @@ rrij=1.0D0/(xj*xj+yj*yj+zj*zj) rij=dsqrt(rrij) sss=sscale(1.0d0/(rij*sigmaii(itypi,itypj))) - sss_ele_cut=sscale_ele(1.0d0/(rij*sigma(itypi,itypj))) - sss_ele_grad=sscagrad_ele(1.0d0/(rij*sigma(itypi,itypj))) + sss_ele_cut=sscale_ele(1.0d0/(rij)) + sss_ele_grad=sscagrad_ele(1.0d0/(rij)) sss_grad=sscale_grad(1.0d0/(rij*sigmaii(itypi,itypj))) if (sss_ele_cut.le.0.0) cycle if (sss.lt.1.0d0) then @@ -14077,7 +14079,7 @@ sigder=fac*sigder fac=rij*fac fac=fac+evdwij*(sss_ele_grad/sss_ele_cut& - /sigma(itypi,itypj)*rij-sss_grad/(1.0-sss)*rij & + *rij-sss_grad/(1.0-sss)*rij & /sigmaii(itypi,itypj)) ! fac=0.0d0 ! Calculate the radial part of the gradient @@ -14311,8 +14313,8 @@ rij=dsqrt(rrij) sss=sscale(1.0d0/(rij*sigmaii(itypi,itypj))) sss_grad=sscale_grad(1.0d0/(rij*sigmaii(itypi,itypj))) - sss_ele_cut=sscale_ele(1.0d0/(rij*sigma(itypi,itypj))) - sss_ele_grad=sscagrad_ele(1.0d0/(rij*sigma(itypi,itypj))) + sss_ele_cut=sscale_ele(1.0d0/(rij)) + sss_ele_grad=sscagrad_ele(1.0d0/(rij)) if (sss_ele_cut.le.0.0) cycle if (sss.gt.0.0d0) then @@ -14368,7 +14370,7 @@ sigder=fac*sigder fac=rij*fac fac=fac+evdwij*(sss_ele_grad/sss_ele_cut& - /sigma(itypi,itypj)*rij+sss_grad/sss*rij & + *rij+sss_grad/sss*rij & /sigmaii(itypi,itypj)) ! fac=0.0d0 diff --git a/source/wham/CMakeLists.txt b/source/wham/CMakeLists.txt index dc01974..a8c8b50 100644 --- a/source/wham/CMakeLists.txt +++ b/source/wham/CMakeLists.txt @@ -86,6 +86,9 @@ if(UNRES_MD_FF STREQUAL "GAB" ) elseif(UNRES_MD_FF STREQUAL "E0LL2Y") # set preprocesor flags set(CPPFLAGS "PROCOR -DSPLITELE -DSCCORPDB" ) +elseif(UNRES_MD_FF STREQUAL "NEWCORR") + # set preprocesor flags + set(CPPFLAGS "PROCOR -DUNRES -DISNAN -DSPLITELE -DLANG0 -DNEWCORR -DCORRCD" ) elseif(UNRES_MD_FF STREQUAL "4P") set(CPPFLAGS "SPLITELE -DLANG0 -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORPDB" ) endif(UNRES_MD_FF STREQUAL "GAB") diff --git a/source/wham/conform_compar.F90 b/source/wham/conform_compar.F90 index 8476646..88c6df9 100644 --- a/source/wham/conform_compar.F90 +++ b/source/wham/conform_compar.F90 @@ -997,9 +997,10 @@ mnum=molnum(i) i1=icont(1,i) i2=icont(2,i) - it1=itype(i1,molnum(i1)) - it2=itype(i2,molnum(i2)) -! print *,"CONTACT",i1,mnum,it1,it2 + it1=itype(i1,1) + it2=itype(i2,1) + if (mnum.eq.0) mnum=1 + print *,"CONTACT",i1,i2,mnum,it1,it2 write (iout,'(i3,2x,a,i4,2x,a,i4,5f8.3,3f10.5)') & i,restyp(it1,mnum),i1,restyp(it2,mnum),i2,cscore(i),& sc_cutoff(iabs(it1),iabs(it2)),ddsc(i),ddla(i),ddlb(i),& @@ -1058,6 +1059,7 @@ i2=icont(2,i) it1=itype(i1,1) it2=itype(i2,1) + write(iout,*) "test",i1,i2,it1,it2 write (iout,'(i3,2x,a,i4,2x,a,i4)') & i,restyp(it1,mnum),i1,restyp(it2,mnum),i2 enddo diff --git a/source/wham/enecalc.F90 b/source/wham/enecalc.F90 index b8fd48f..b175378 100644 --- a/source/wham/enecalc.F90 +++ b/source/wham/enecalc.F90 @@ -2,6 +2,7 @@ !----------------------------------------------------------------------------- use io_units use wham_data + use control_data, only: tor_mode ! use geometry_data, only:nres,boxxsize,boxysize,boxzsize use energy_data @@ -675,6 +676,8 @@ write(iout,*)"enecalc_ i ntot",i,ntot ww_all(16,iparm)=wvdwpp ww_all(17,iparm)=wbond ww_all(19,iparm)=wsccor + ww_all(42,iparm)=wcatprot + ww_all(41,iparm)=wcatcat ! Store bond parameters vbldp0_all(iparm)=vbldp0 akp_all(iparm)=akp @@ -709,6 +712,7 @@ write(iout,*)"enecalc_ i ntot",i,ntot sigc0_all(i,iparm)=sigc0(i) enddo #else + IF (tor_mode.eq.0) THEN nthetyp_all(iparm)=nthetyp ntheterm_all(iparm)=ntheterm ntheterm2_all(iparm)=ntheterm2 @@ -760,6 +764,9 @@ write(iout,*)"enecalc_ i ntot",i,ntot enddo enddo enddo + ELSE + write(iout,*) "Need storing parameters" + ENDIF #endif #ifdef CRYST_SC ! Store the sidechain rotamer parameters @@ -787,6 +794,7 @@ write(iout,*)"enecalc_ i ntot",i,ntot enddo enddo #endif + IF (tor_mode.eq.0) THEN ! Store the torsional parameters do iblock=1,2 do i=-ntortyp+1,ntortyp-1 @@ -845,6 +853,9 @@ write(iout,*)"enecalc_ i ntot",i,ntot enddo enddo enddo + ELSE + write(iout,*) "NEED storing parameters" + ENDIF ! Store the parameters of electrostatic interactions do i=1,2 do j=1,2 @@ -947,6 +958,8 @@ write(iout,*)"end of store_parm" wvdwpp=ww_all(16,iparm) wbond=ww_all(17,iparm) wsccor=ww_all(19,iparm) + wcatprot=ww_all(42,iparm) + wcatcat=ww_all(41,iparm) ! Restore bond parameters vbldp0=vbldp0_all(iparm) akp=akp_all(iparm) @@ -981,6 +994,7 @@ write(iout,*)"end of store_parm" sigc0(i)=sigc0_all(i,iparm) enddo #else + IF (tor_mode.eq.0) THEN nthetyp=nthetyp_all(iparm) ntheterm=ntheterm_all(iparm) ntheterm2=ntheterm2_all(iparm) @@ -1032,6 +1046,9 @@ write(iout,*)"end of store_parm" enddo enddo enddo + ELSE + write (iout,*) "Need storing parameters" + ENDIF #endif ! Restore the sidechain rotamer parameters #ifdef CRYST_SC @@ -1058,6 +1075,7 @@ write(iout,*)"end of store_parm" enddo enddo #endif + IF (tor_mode.eq.0) THEN ! Restore the torsional parameters do iblock=1,2 do i=-ntortyp+1,ntortyp-1 @@ -1099,6 +1117,7 @@ write(iout,*)"end of store_parm" enddo enddo enddo + ! Restore parameters of the cumulants do i=-nloctyp,nloctyp do j=1,2 @@ -1116,6 +1135,9 @@ write(iout,*)"end of store_parm" enddo enddo enddo + ELSE + write(iout,*) "need storing parameters" + ENDIF ! Restore the parameters of electrostatic interactions do i=1,2 do j=1,2 diff --git a/source/wham/io_wham.F90 b/source/wham/io_wham.F90 index 7ebde6f..958e244 100644 --- a/source/wham/io_wham.F90 +++ b/source/wham/io_wham.F90 @@ -332,7 +332,7 @@ use geometry_data use energy_data use control_data, only: maxterm,maxlor,maxterm_sccor,& !maxtor - maxtermd_1,maxtermd_2 !,maxthetyp,maxthetyp1 + maxtermd_1,maxtermd_2,tor_mode,scelemode !,maxthetyp,maxthetyp1 use MD_data !el use MPI_data !el use map_data @@ -365,7 +365,7 @@ character(len=1) :: t1,t2,t3 character(len=1) :: onelett(4) = (/"G","A","P","D"/) character(len=1) :: toronelet(-2:2) = (/"p","a","G","A","P"/) - logical :: lprint + logical :: lprint,SPLIT_FOURIERTOR real(kind=8),dimension(3,3,maxlob) :: blower !(3,3,maxlob) character(len=800) :: controlcard character(len=256) :: bondname_t,thetname_t,rotname_t,torname_t,& @@ -374,14 +374,16 @@ !el integer ilen !el external ilen character(len=16) :: key - integer :: iparm + integer :: iparm,nkcctyp !el real(kind=8) :: ip,mp real(kind=8) :: dwa16,akl,si,rri,epsij,rrij,sigeps,sigt1sq,epsijlip,& sigt2sq,sigii1,sigii2,ratsig1,ratsig2,rsum_max,r_augm real(kind=8) :: v0ij,v0ijsccor,v0ijsccor1,v0ijsccor2,v0ijsccor3,rjunk,& res1 - integer :: i,j,ichir1,ichir2,k,l,m,kk,ii,mm,junk,lll,ll,llll,n + integer :: i,j,ichir1,ichir2,k,l,m,kk,ii,mm,junk,lll,ll,llll,n,jj integer :: nlobi,iblock,maxinter,iscprol,ncatprotparm + character*3 string + ! ! Body ! @@ -428,8 +430,8 @@ allocate(ww(max_eneW)) wvdwpp=ww(16) wbond=ww(18) wsccor=ww(19) - wcatcat=ww(44) - wcatprot=ww(43) + wcatcat=ww(42) + wcatprot=ww(41) endif ! @@ -454,8 +456,8 @@ allocate(ww(max_eneW)) weights(17)=wbond weights(18)=0 !scal14 ! weights(21)=wsccor - weights(43)=wcatprot - weights(44)=wcatcat + weights(42)=wcatprot + weights(41)=wcatcat ! el-------- call card_concat(controlcard,.false.) @@ -738,13 +740,16 @@ allocate(ww(max_eneW)) ! Read the parameters of Utheta determined from ab initio surfaces ! Kozlowska et al., J. Phys.: Condens. Matter 19 (2007) 285203 ! + allocate(ithetyp(-ntyp1:ntyp1)) !(-ntyp1:ntyp1) ! write (iout,*) "tu dochodze" + IF (tor_mode.eq.0) THEN +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! read (ithep,*) nthetyp,ntheterm,ntheterm2,& ntheterm3,nsingle,ndouble nntheterm=max0(ntheterm,ntheterm2,ntheterm3) !---------------------------------------------------- - allocate(ithetyp(-ntyp1:ntyp1)) !(-ntyp1:ntyp1) +! allocate(ithetyp(-ntyp1:ntyp1)) !(-ntyp1:ntyp1) allocate(aa0thet(-nthetyp-1:nthetyp+1,& -nthetyp-1:nthetyp+1,-nthetyp-1:nthetyp+1,2)) !(-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2) @@ -864,7 +869,7 @@ allocate(ww(max_eneW)) ! ! Control printout of the coefficients of virtual-bond-angle potentials ! -do iblock=1,2 + do iblock=1,2 if (lprint) then write (iout,'(//a)') 'Parameter of virtual-bond-angle potential' do i=1,nthetyp+1 @@ -903,7 +908,30 @@ do iblock=1,2 enddo call flush(iout) endif -enddo + enddo + ELSE +!C here will be the apropriate recalibrating for D-aminoacid + read (ithep,*) nthetyp + allocate(nbend_kcc_Tb(-nthetyp:nthetyp)) + allocate(v1bend_chyb(0:36,-nthetyp:nthetyp)) + do i=-nthetyp+1,nthetyp-1 + read (ithep,*) nbend_kcc_Tb(i) + do j=0,nbend_kcc_Tb(i) + read (ithep,*) ijunk,v1bend_chyb(j,i) + enddo + enddo + if (lprint) then + write (iout,'(a)') & + "Parameters of the valence-only potentials" + do i=-nthetyp+1,nthetyp-1 + write (iout,'(2a)') "Type ",toronelet(i) + do k=0,nbend_kcc_Tb(i) + write(iout,'(i5,f15.5)') k,v1bend_chyb(k,i) + enddo + enddo + endif + ENDIF ! TOR_MODE + #endif !------------------------------------------- allocate(nlob(ntyp1)) !(ntyp1) @@ -1042,6 +1070,358 @@ enddo enddo #endif close(irotam) + + read (ifourier,*) nloctyp +!el write(iout,*)"nloctyp",nloctyp + SPLIT_FOURIERTOR = nloctyp.lt.0 + nloctyp = iabs(nloctyp) +#ifdef NEWCORR + if (.not.allocated(itype2loc)) allocate(itype2loc(-ntyp1:ntyp1)) + print *,"shape",shape(itype2loc) + allocate(iloctyp(-nloctyp:nloctyp)) + allocate(bnew1(3,2,-nloctyp:nloctyp)) + allocate(bnew2(3,2,-nloctyp:nloctyp)) + allocate(ccnew(3,2,-nloctyp:nloctyp)) + allocate(ddnew(3,2,-nloctyp:nloctyp)) + allocate(e0new(3,-nloctyp:nloctyp)) + allocate(eenew(2,2,2,-nloctyp:nloctyp)) + allocate(bnew1tor(3,2,-nloctyp:nloctyp)) + allocate(bnew2tor(3,2,-nloctyp:nloctyp)) + allocate(ccnewtor(3,2,-nloctyp:nloctyp)) + allocate(ddnewtor(3,2,-nloctyp:nloctyp)) + allocate(e0newtor(3,-nloctyp:nloctyp)) + allocate(eenewtor(2,2,2,-nloctyp:nloctyp)) + + read (ifourier,*) (itype2loc(i),i=1,ntyp) + read (ifourier,*) (iloctyp(i),i=0,nloctyp-1) + itype2loc(ntyp1)=nloctyp + iloctyp(nloctyp)=ntyp1 + do i=1,ntyp1 + itype2loc(-i)=-itype2loc(i) + enddo +#else + allocate(iloctyp(-nloctyp:nloctyp)) + allocate(itype2loc(-ntyp1:ntyp1)) + iloctyp(0)=10 + iloctyp(1)=9 + iloctyp(2)=20 + iloctyp(3)=ntyp1 +#endif + do i=1,nloctyp + iloctyp(-i)=-iloctyp(i) + enddo +!c write (iout,*) "itype2loc",(itype2loc(i),i=1,ntyp1) +!c write (iout,*) "nloctyp",nloctyp, +!c & " iloctyp",(iloctyp(i),i=0,nloctyp) +!c write (iout,*) "itype2loc",(itype2loc(i),i=1,ntyp1) +!c write (iout,*) "nloctyp",nloctyp, +!c & " iloctyp",(iloctyp(i),i=0,nloctyp) +#ifdef NEWCORR + do i=0,nloctyp-1 +!c write (iout,*) "NEWCORR",i + read (ifourier,*) + do ii=1,3 + do j=1,2 + read (ifourier,*) bnew1(ii,j,i) + enddo + enddo +!c write (iout,*) "NEWCORR BNEW1" +!c write (iout,*) ((bnew1(ii,j,i),ii=1,3),j=1,2) + do ii=1,3 + do j=1,2 + read (ifourier,*) bnew2(ii,j,i) + enddo + enddo +!c write (iout,*) "NEWCORR BNEW2" +!c write (iout,*) ((bnew2(ii,j,i),ii=1,3),j=1,2) + do kk=1,3 + read (ifourier,*) ccnew(kk,1,i) + read (ifourier,*) ccnew(kk,2,i) + enddo +!c write (iout,*) "NEWCORR CCNEW" +!c write (iout,*) ((ccnew(ii,j,i),ii=1,3),j=1,2) + do kk=1,3 + read (ifourier,*) ddnew(kk,1,i) + read (ifourier,*) ddnew(kk,2,i) + enddo +!c write (iout,*) "NEWCORR DDNEW" +!c write (iout,*) ((ddnew(ii,j,i),ii=1,3),j=1,2) + do ii=1,2 + do jj=1,2 + do kk=1,2 + read (ifourier,*) eenew(ii,jj,kk,i) + enddo + enddo + enddo +!c write (iout,*) "NEWCORR EENEW1" +!c write(iout,*)(((eenew(ii,jj,kk,i),kk=1,2),jj=1,2),ii=1,2) + do ii=1,3 + read (ifourier,*) e0new(ii,i) + enddo +!c write (iout,*) (e0new(ii,i),ii=1,3) + enddo +!c write (iout,*) "NEWCORR EENEW" + do i=0,nloctyp-1 + do ii=1,3 + ccnew(ii,1,i)=ccnew(ii,1,i)/2 + ccnew(ii,2,i)=ccnew(ii,2,i)/2 + ddnew(ii,1,i)=ddnew(ii,1,i)/2 + ddnew(ii,2,i)=ddnew(ii,2,i)/2 + enddo + enddo + do i=1,nloctyp-1 + do ii=1,3 + bnew1(ii,1,-i)= bnew1(ii,1,i) + bnew1(ii,2,-i)=-bnew1(ii,2,i) + bnew2(ii,1,-i)= bnew2(ii,1,i) + bnew2(ii,2,-i)=-bnew2(ii,2,i) + enddo + do ii=1,3 +!c ccnew(ii,1,i)=ccnew(ii,1,i)/2 +!c ccnew(ii,2,i)=ccnew(ii,2,i)/2 +!c ddnew(ii,1,i)=ddnew(ii,1,i)/2 +!c ddnew(ii,2,i)=ddnew(ii,2,i)/2 + ccnew(ii,1,-i)=ccnew(ii,1,i) + ccnew(ii,2,-i)=-ccnew(ii,2,i) + ddnew(ii,1,-i)=ddnew(ii,1,i) + ddnew(ii,2,-i)=-ddnew(ii,2,i) + enddo + e0new(1,-i)= e0new(1,i) + e0new(2,-i)=-e0new(2,i) + e0new(3,-i)=-e0new(3,i) + do kk=1,2 + eenew(kk,1,1,-i)= eenew(kk,1,1,i) + eenew(kk,1,2,-i)=-eenew(kk,1,2,i) + eenew(kk,2,1,-i)=-eenew(kk,2,1,i) + eenew(kk,2,2,-i)= eenew(kk,2,2,i) + enddo + enddo + if (lprint) then + write (iout,'(a)') "Coefficients of the multibody terms" + do i=-nloctyp+1,nloctyp-1 + write (iout,*) "Type: ",onelet(iloctyp(i)) + write (iout,*) "Coefficients of the expansion of B1" + do j=1,2 + write (iout,'(3hB1(,i1,1h),3f10.5)') j,(bnew1(k,j,i),k=1,3) + enddo + write (iout,*) "Coefficients of the expansion of B2" + do j=1,2 + write (iout,'(3hB2(,i1,1h),3f10.5)') j,(bnew2(k,j,i),k=1,3) + enddo + write (iout,*) "Coefficients of the expansion of C" + write (iout,'(3hC11,3f10.5)') (ccnew(j,1,i),j=1,3) + write (iout,'(3hC12,3f10.5)') (ccnew(j,2,i),j=1,3) + write (iout,*) "Coefficients of the expansion of D" + write (iout,'(3hD11,3f10.5)') (ddnew(j,1,i),j=1,3) + write (iout,'(3hD12,3f10.5)') (ddnew(j,2,i),j=1,3) + write (iout,*) "Coefficients of the expansion of E" + write (iout,'(2hE0,3f10.5)') (e0new(j,i),j=1,3) + do j=1,2 + do k=1,2 + write (iout,'(1hE,2i1,2f10.5)') j,k,(eenew(l,j,k,i),l=1,2) + enddo + enddo + enddo + endif + IF (SPLIT_FOURIERTOR) THEN + do i=0,nloctyp-1 +!c write (iout,*) "NEWCORR TOR",i + read (ifourier,*) + do ii=1,3 + do j=1,2 + read (ifourier,*) bnew1tor(ii,j,i) + enddo + enddo +!c write (iout,*) "NEWCORR BNEW1 TOR" +!c write (iout,*) ((bnew1tor(ii,j,i),ii=1,3),j=1,2) + do ii=1,3 + do j=1,2 + read (ifourier,*) bnew2tor(ii,j,i) + enddo + enddo +!c write (iout,*) "NEWCORR BNEW2 TOR" +!c write (iout,*) ((bnew2tor(ii,j,i),ii=1,3),j=1,2) + do kk=1,3 + read (ifourier,*) ccnewtor(kk,1,i) + read (ifourier,*) ccnewtor(kk,2,i) + enddo +!c write (iout,*) "NEWCORR CCNEW TOR" +!c write (iout,*) ((ccnew(ii,j,i),ii=1,3),j=1,2) + do kk=1,3 + read (ifourier,*) ddnewtor(kk,1,i) + read (ifourier,*) ddnewtor(kk,2,i) + enddo +!c write (iout,*) "NEWCORR DDNEW TOR" +!c write (iout,*) ((ddnewtor(ii,j,i),ii=1,3),j=1,2) + do ii=1,2 + do jj=1,2 + do kk=1,2 + read (ifourier,*) eenewtor(ii,jj,kk,i) + enddo + enddo + enddo +!c write (iout,*) "NEWCORR EENEW1 TOR" +!c write(iout,*)(((eenewtor(ii,jj,kk,i),kk=1,2),jj=1,2),ii=1,2) + do ii=1,3 + read (ifourier,*) e0newtor(ii,i) + enddo +!c write (iout,*) (e0newtor(ii,i),ii=1,3) + enddo +!c write (iout,*) "NEWCORR EENEW TOR" + do i=0,nloctyp-1 + do ii=1,3 + ccnewtor(ii,1,i)=ccnewtor(ii,1,i)/2 + ccnewtor(ii,2,i)=ccnewtor(ii,2,i)/2 + ddnewtor(ii,1,i)=ddnewtor(ii,1,i)/2 + ddnewtor(ii,2,i)=ddnewtor(ii,2,i)/2 + enddo + enddo + do i=1,nloctyp-1 + do ii=1,3 + bnew1tor(ii,1,-i)= bnew1tor(ii,1,i) + bnew1tor(ii,2,-i)=-bnew1tor(ii,2,i) + bnew2tor(ii,1,-i)= bnew2tor(ii,1,i) + bnew2tor(ii,2,-i)=-bnew2tor(ii,2,i) + enddo + do ii=1,3 + ccnewtor(ii,1,-i)=ccnewtor(ii,1,i) + ccnewtor(ii,2,-i)=-ccnewtor(ii,2,i) + ddnewtor(ii,1,-i)=ddnewtor(ii,1,i) + ddnewtor(ii,2,-i)=-ddnewtor(ii,2,i) + enddo + e0newtor(1,-i)= e0newtor(1,i) + e0newtor(2,-i)=-e0newtor(2,i) + e0newtor(3,-i)=-e0newtor(3,i) + do kk=1,2 + eenewtor(kk,1,1,-i)= eenewtor(kk,1,1,i) + eenewtor(kk,1,2,-i)=-eenewtor(kk,1,2,i) + eenewtor(kk,2,1,-i)=-eenewtor(kk,2,1,i) + eenewtor(kk,2,2,-i)= eenewtor(kk,2,2,i) + enddo + enddo + if (lprint) then + write (iout,'(a)') & + "Single-body coefficients of the torsional potentials" + do i=-nloctyp+1,nloctyp-1 + write (iout,*) "Type: ",onelet(iloctyp(i)) + write (iout,*) "Coefficients of the expansion of B1tor" + do j=1,2 + write (iout,'(3hB1(,i1,1h),3f10.5)') & + j,(bnew1tor(k,j,i),k=1,3) + enddo + write (iout,*) "Coefficients of the expansion of B2tor" + do j=1,2 + write (iout,'(3hB2(,i1,1h),3f10.5)') & + j,(bnew2tor(k,j,i),k=1,3) + enddo + write (iout,*) "Coefficients of the expansion of Ctor" + write (iout,'(3hC11,3f10.5)') (ccnewtor(j,1,i),j=1,3) + write (iout,'(3hC12,3f10.5)') (ccnewtor(j,2,i),j=1,3) + write (iout,*) "Coefficients of the expansion of Dtor" + write (iout,'(3hD11,3f10.5)') (ddnewtor(j,1,i),j=1,3) + write (iout,'(3hD12,3f10.5)') (ddnewtor(j,2,i),j=1,3) + write (iout,*) "Coefficients of the expansion of Etor" + write (iout,'(2hE0,3f10.5)') (e0newtor(j,i),j=1,3) + do j=1,2 + do k=1,2 + write (iout,'(1hE,2i1,2f10.5)') & + j,k,(eenewtor(l,j,k,i),l=1,2) + enddo + enddo + enddo + endif + ELSE + do i=-nloctyp+1,nloctyp-1 + do ii=1,3 + do j=1,2 + bnew1tor(ii,j,i)=bnew1(ii,j,i) + enddo + enddo + do ii=1,3 + do j=1,2 + bnew2tor(ii,j,i)=bnew2(ii,j,i) + enddo + enddo + do ii=1,3 + ccnewtor(ii,1,i)=ccnew(ii,1,i) + ccnewtor(ii,2,i)=ccnew(ii,2,i) + ddnewtor(ii,1,i)=ddnew(ii,1,i) + ddnewtor(ii,2,i)=ddnew(ii,2,i) + enddo + enddo + ENDIF !SPLIT_FOURIER_TOR +#else + allocate(ccold(2,2,-nloctyp-1:nloctyp+1)) + allocate(ddold(2,2,-nloctyp-1:nloctyp+1)) + allocate(eeold(2,2,-nloctyp-1:nloctyp+1)) + allocate(b(13,-nloctyp-1:nloctyp+1)) + if (lprint) & + write (iout,*) "Coefficients of the expansion of Eloc(l1,l2)" + do i=0,nloctyp-1 + read (ifourier,*) + read (ifourier,*) (b(ii,i),ii=1,13) + if (lprint) then + write (iout,*) 'Type ',onelet(iloctyp(i)) + write (iout,'(a,i2,a,f10.5)') ('b(',ii,')=',b(ii,i),ii=1,13) + endif + if (i.gt.0) then + b(2,-i)= b(2,i) + b(3,-i)= b(3,i) + b(4,-i)=-b(4,i) + b(5,-i)=-b(5,i) + endif + CCold(1,1,i)= b(7,i) + CCold(2,2,i)=-b(7,i) + CCold(2,1,i)= b(9,i) + CCold(1,2,i)= b(9,i) + CCold(1,1,-i)= b(7,i) + CCold(2,2,-i)=-b(7,i) + CCold(2,1,-i)=-b(9,i) + CCold(1,2,-i)=-b(9,i) + DDold(1,1,i)= b(6,i) + DDold(2,2,i)=-b(6,i) + DDold(2,1,i)= b(8,i) + DDold(1,2,i)= b(8,i) + DDold(1,1,-i)= b(6,i) + DDold(2,2,-i)=-b(6,i) + DDold(2,1,-i)=-b(8,i) + DDold(1,2,-i)=-b(8,i) + EEold(1,1,i)= b(10,i)+b(11,i) + EEold(2,2,i)=-b(10,i)+b(11,i) + EEold(2,1,i)= b(12,i)-b(13,i) + EEold(1,2,i)= b(12,i)+b(13,i) + EEold(1,1,-i)= b(10,i)+b(11,i) + EEold(2,2,-i)=-b(10,i)+b(11,i) + EEold(2,1,-i)=-b(12,i)+b(13,i) + EEold(1,2,-i)=-b(12,i)-b(13,i) + write(iout,*) "TU DOCHODZE" + print *,"JESTEM" + enddo + if (lprint) then + write (iout,*) + write (iout,*) & + "Coefficients of the cumulants (independent of valence angles)" + do i=-nloctyp+1,nloctyp-1 + write (iout,*) 'Type ',onelet(iloctyp(i)) + write (iout,*) 'B1' + write(iout,'(2f10.5)') B(3,i),B(5,i) + write (iout,*) 'B2' + write(iout,'(2f10.5)') B(2,i),B(4,i) + write (iout,*) 'CC' + do j=1,2 + write (iout,'(2f10.5)') CCold(j,1,i),CCold(j,2,i) + enddo + write(iout,*) 'DD' + do j=1,2 + write (iout,'(2f10.5)') DDold(j,1,i),DDold(j,2,i) + enddo + write(iout,*) 'EE' + do j=1,2 + write (iout,'(2f10.5)') EEold(j,1,i),EEold(j,2,i) + enddo + enddo + endif +#endif #ifdef CRYST_TOR ! ! Read torsional parameters in old format @@ -1081,6 +1461,8 @@ enddo ! ! Read torsional parameters ! + IF (TOR_MODE.eq.0) THEN + allocate(itortyp(-ntyp1:ntyp1)) !(-ntyp1:ntyp1) read (itorp,*) ntortyp @@ -1269,6 +1651,113 @@ enddo enddo enddo endif + ELSE IF (TOR_MODE.eq.1) THEN + +!C read valence-torsional parameters + read (itorp,*) ntortyp + nkcctyp=ntortyp + write (iout,*) "Valence-torsional parameters read in ntortyp",& + ntortyp + read (itorp,*) (itortyp(i),i=1,ntyp) + write (iout,*) "itortyp_kcc",(itortyp(i),i=1,ntyp) +#ifndef NEWCORR + do i=1,ntyp1 + itype2loc(i)=itortyp(i) + enddo +#endif + do i=-ntyp,-1 + itortyp(i)=-itortyp(-i) + enddo + do i=-ntortyp+1,ntortyp-1 + do j=-ntortyp+1,ntortyp-1 +!C first we read the cos and sin gamma parameters + read (itorp,'(13x,a)') string + write (iout,*) i,j,string + read (itorp,*) & + nterm_kcc(j,i),nterm_kcc_Tb(j,i) +!C read (itorkcc,*,end=121,err=121) nterm_kcc_Tb(j,i) + do k=1,nterm_kcc(j,i) + do l=1,nterm_kcc_Tb(j,i) + do ll=1,nterm_kcc_Tb(j,i) + read (itorp,*) ii,jj,kk, & + v1_kcc(ll,l,k,j,i),v2_kcc(ll,l,k,j,i) + enddo + enddo + enddo + enddo + enddo + ELSE +#ifdef NEWCORR +!c AL 4/8/16: Calculate coefficients from one-body parameters + ntortyp=nloctyp + allocate(itortyp(-ntyp1:ntyp1)) !(-ntyp1:ntyp1) + allocate(nterm_kcc(-ntyp1:ntyp1,-ntyp1:ntyp1)) + allocate(nterm_kcc_Tb(-ntyp1:ntyp1,-ntyp1:ntyp1)) + allocate(v1_kcc(6,6,6,-ntyp1:ntyp1,-ntyp1:ntyp1)) + allocate(v2_kcc(6,6,6,-ntyp1:ntyp1,-ntyp1:ntyp1)) + + do i=-ntyp1,ntyp1 + print *,i,itortyp(i) + itortyp(i)=itype2loc(i) + enddo + write (iout,*) & + "Val-tor parameters calculated from cumulant coefficients ntortyp"& + ,ntortyp + do i=-ntortyp+1,ntortyp-1 + do j=-ntortyp+1,ntortyp-1 + nterm_kcc(j,i)=2 + nterm_kcc_Tb(j,i)=3 + do k=1,nterm_kcc_Tb(j,i) + do l=1,nterm_kcc_Tb(j,i) + v1_kcc(k,l,1,i,j)=bnew1tor(k,1,i)*bnew2tor(l,1,j)& + +bnew1tor(k,2,i)*bnew2tor(l,2,j) + v2_kcc(k,l,1,i,j)=bnew1tor(k,1,i)*bnew2tor(l,2,j)& + +bnew1tor(k,2,i)*bnew2tor(l,1,j) + enddo + enddo + do k=1,nterm_kcc_Tb(j,i) + do l=1,nterm_kcc_Tb(j,i) +#ifdef CORRCD + v1_kcc(k,l,2,i,j)=-(ccnewtor(k,1,i)*ddnewtor(l,1,j) & + -ccnewtor(k,2,i)*ddnewtor(l,2,j)) + v2_kcc(k,l,2,i,j)=-(ccnewtor(k,2,i)*ddnewtor(l,1,j) & + +ccnewtor(k,1,i)*ddnewtor(l,2,j)) +#else + v1_kcc(k,l,2,i,j)=-0.25*(ccnewtor(k,1,i)*ddnewtor(l,1,j) & + -ccnewtor(k,2,i)*ddnewtor(l,2,j)) + v2_kcc(k,l,2,i,j)=-0.25*(ccnewtor(k,2,i)*ddnewtor(l,1,j) & + +ccnewtor(k,1,i)*ddnewtor(l,2,j)) +#endif + enddo + enddo +!c f(theta,gamma)=-(b21(theta)*b11(theta)+b12(theta)*b22(theta))*cos(gamma)-(b22(theta)*b11(theta)+b21(theta)*b12(theta))*sin(gamma)+(c11(theta)*d11(theta)-c12(theta)*d12(theta))*cos(2*gamma)+(c12(theta)*d11(theta)+c11(theta)*d12(theta))*sin(2*gamma) + enddo + enddo +#else + write (iout,*) "TOR_MODE>1 only with NEWCORR" + stop +#endif + ENDIF ! TOR_MODE + if (tor_mode.gt.0 .and. lprint) then +!c Print valence-torsional parameters + write (iout,'(a)') & + "Parameters of the valence-torsional potentials" + do i=-ntortyp+1,ntortyp-1 + do j=-ntortyp+1,ntortyp-1 + write (iout,'(3a)') "Type ",toronelet(i),toronelet(j) + write (iout,'(3a5,2a15)') "itor","ival","jval","v_kcc","v2_kcc" + do k=1,nterm_kcc(j,i) + do l=1,nterm_kcc_Tb(j,i) + do ll=1,nterm_kcc_Tb(j,i) + write (iout,'(3i5,2f15.4)')& + k,l-1,ll-1,v1_kcc(ll,l,k,j,i),v2_kcc(ll,l,k,j,i) + enddo + enddo + enddo + enddo + enddo + endif + #endif !elwrite(iout,*) "parmread kontrol sc-bb" ! Read of Side-chain backbone correlation parameters @@ -1397,153 +1886,7 @@ enddo enddo enddo endif -! -! 9/18/99 (AL) Read coefficients of the Fourier expansion of the local -! interaction energy of the Gly, Ala, and Pro prototypes. -! - read (ifourier,*) nloctyp -!el write(iout,*)"nloctyp",nloctyp -!el from module energy------- - allocate(b1(2,-nloctyp-1:nloctyp+1)) !(2,-maxtor:maxtor) - allocate(b2(2,-nloctyp-1:nloctyp+1)) !(2,-maxtor:maxtor) - allocate(b1tilde(2,-nloctyp-1:nloctyp+1)) !(2,-maxtor:maxtor) - allocate(cc(2,2,-nloctyp-1:nloctyp+1)) - allocate(dd(2,2,-nloctyp-1:nloctyp+1)) - allocate(ee(2,2,-nloctyp-1:nloctyp+1)) - allocate(ctilde(2,2,-nloctyp-1:nloctyp+1)) - allocate(dtilde(2,2,-nloctyp-1:nloctyp+1)) !(2,2,-maxtor:maxtor) - do i=1,2 - do ii=-nloctyp-1,nloctyp+1 - b1(i,ii)=0.0d0 - b2(i,ii)=0.0d0 - b1tilde(i,ii)=0.0d0 - do j=1,2 - cc(j,i,ii)=0.0d0 - dd(j,i,ii)=0.0d0 - ee(j,i,ii)=0.0d0 - ctilde(j,i,ii)=0.0d0 - dtilde(j,i,ii)=0.0d0 - enddo - enddo - enddo -!-------------------------------- - allocate(b(13,0:nloctyp)) - - do i=0,nloctyp-1 - read (ifourier,*) - read (ifourier,*) (b(ii,i),ii=1,13) - if (lprint) then - write (iout,*) 'Type',i - write (iout,'(a,i2,a,f10.5)') ('b(',ii,')=',b(ii,i),ii=1,13) - endif - B1(1,i) = b(3,i) - B1(2,i) = b(5,i) - B1(1,-i) = b(3,i) - B1(2,-i) = -b(5,i) -! b1(1,i)=0.0d0 -! b1(2,i)=0.0d0 - B1tilde(1,i) = b(3,i) - B1tilde(2,i) =-b(5,i) - B1tilde(1,-i) =-b(3,i) - B1tilde(2,-i) =b(5,i) -! b1tilde(1,i)=0.0d0 -! b1tilde(2,i)=0.0d0 - B2(1,i) = b(2,i) - B2(2,i) = b(4,i) - B2(1,-i) =b(2,i) - B2(2,-i) =-b(4,i) - -! b2(1,i)=0.0d0 -! b2(2,i)=0.0d0 - CC(1,1,i)= b(7,i) - CC(2,2,i)=-b(7,i) - CC(2,1,i)= b(9,i) - CC(1,2,i)= b(9,i) - CC(1,1,-i)= b(7,i) - CC(2,2,-i)=-b(7,i) - CC(2,1,-i)=-b(9,i) - CC(1,2,-i)=-b(9,i) -! CC(1,1,i)=0.0d0 -! CC(2,2,i)=0.0d0 -! CC(2,1,i)=0.0d0 -! CC(1,2,i)=0.0d0 - Ctilde(1,1,i)=b(7,i) - Ctilde(1,2,i)=b(9,i) - Ctilde(2,1,i)=-b(9,i) - Ctilde(2,2,i)=b(7,i) - Ctilde(1,1,-i)=b(7,i) - Ctilde(1,2,-i)=-b(9,i) - Ctilde(2,1,-i)=b(9,i) - Ctilde(2,2,-i)=b(7,i) - -! Ctilde(1,1,i)=0.0d0 -! Ctilde(1,2,i)=0.0d0 -! Ctilde(2,1,i)=0.0d0 -! Ctilde(2,2,i)=0.0d0 - DD(1,1,i)= b(6,i) - DD(2,2,i)=-b(6,i) - DD(2,1,i)= b(8,i) - DD(1,2,i)= b(8,i) - DD(1,1,-i)= b(6,i) - DD(2,2,-i)=-b(6,i) - DD(2,1,-i)=-b(8,i) - DD(1,2,-i)=-b(8,i) -! DD(1,1,i)=0.0d0 -! DD(2,2,i)=0.0d0 -! DD(2,1,i)=0.0d0 -! DD(1,2,i)=0.0d0 - Dtilde(1,1,i)=b(6,i) - Dtilde(1,2,i)=b(8,i) - Dtilde(2,1,i)=-b(8,i) - Dtilde(2,2,i)=b(6,i) - Dtilde(1,1,-i)=b(6,i) - Dtilde(1,2,-i)=-b(8,i) - Dtilde(2,1,-i)=b(8,i) - Dtilde(2,2,-i)=b(6,i) - -! Dtilde(1,1,i)=0.0d0 -! Dtilde(1,2,i)=0.0d0 -! Dtilde(2,1,i)=0.0d0 -! Dtilde(2,2,i)=0.0d0 - EE(1,1,i)= b(10,i)+b(11,i) - EE(2,2,i)=-b(10,i)+b(11,i) - EE(2,1,i)= b(12,i)-b(13,i) - EE(1,2,i)= b(12,i)+b(13,i) - EE(1,1,-i)= b(10,i)+b(11,i) - EE(2,2,-i)=-b(10,i)+b(11,i) - EE(2,1,-i)=-b(12,i)+b(13,i) - EE(1,2,-i)=-b(12,i)-b(13,i) - -! ee(1,1,i)=1.0d0 -! ee(2,2,i)=1.0d0 -! ee(2,1,i)=0.0d0 -! ee(1,2,i)=0.0d0 -! ee(2,1,i)=ee(1,2,i) - enddo - if (lprint) then - do i=1,nloctyp - write (iout,*) 'Type',i - write (iout,*) 'B1' -! write (iout,'(f10.5)') B1(:,i) - write(iout,*) B1(1,i),B1(2,i) - write (iout,*) 'B2' -! write (iout,'(f10.5)') B2(:,i) - write(iout,*) B2(1,i),B2(2,i) - write (iout,*) 'CC' - do j=1,2 - write (iout,'(2f10.5)') CC(j,1,i),CC(j,2,i) - enddo - write(iout,*) 'DD' - do j=1,2 - write (iout,'(2f10.5)') DD(j,1,i),DD(j,2,i) - enddo - write(iout,*) 'EE' - do j=1,2 - write (iout,'(2f10.5)') EE(j,1,i),EE(j,2,i) - enddo - enddo - endif ! ! Read electrostatic-interaction parameters ! @@ -1632,6 +1975,7 @@ enddo !---------------------- GB or BP potential ----------------------------- case (3:4) ! 30 do i=1,ntyp + if (scelemode.eq.0) then do i=1,ntyp read (isidep,*)(eps(i,j),j=i,ntyp) enddo @@ -1658,6 +2002,86 @@ enddo write (iout,'(a3,6x,4f10.5)') (restyp(i,molnum(i)),sigma0(i),sigii(i),& chip(i),alp(i),i=1,ntyp) endif + else + allocate(icharge(ntyp1)) +! print *,ntyp,icharge(i) + icharge(:)=0 + read (isidep,*) (icharge(i),i=1,ntyp) + print *,ntyp,icharge(i) +! if(.not.allocated(eps)) allocate(eps(-ntyp + write (2,*) "icharge",(icharge(i),i=1,ntyp) + allocate(alphapol(ntyp,ntyp),epshead(ntyp,ntyp),sig0head(ntyp,ntyp)) + allocate(sigiso1(ntyp,ntyp),rborn(ntyp,ntyp),sigmap1(ntyp,ntyp)) + allocate(sigmap2(ntyp,ntyp),sigiso2(ntyp,ntyp)) + allocate(chis(ntyp,ntyp),wquad(ntyp,ntyp),chipp(ntyp,ntyp)) + allocate(epsintab(ntyp,ntyp)) + allocate(dtail(2,ntyp,ntyp)) + print *,"control line 1" + allocate(alphasur(4,ntyp,ntyp),alphiso(4,ntyp,ntyp)) + allocate(wqdip(2,ntyp,ntyp)) + allocate(wstate(4,ntyp,ntyp)) + allocate(dhead(2,2,ntyp,ntyp)) + allocate(nstate(ntyp,ntyp)) + allocate(debaykap(ntyp,ntyp)) + print *,"control line 2" + if (.not.allocated(sigma)) allocate(sigma(0:ntyp1,0:ntyp1)) + if (.not.allocated(chi)) allocate(chi(ntyp1,ntyp1)) !(ntyp,ntyp) + + do i=1,ntyp + do j=1,i +! write (*,*) "Im in ALAB", i, " ", j + read(isidep,*) & + eps(i,j),sigma(i,j),chi(i,j),chi(j,i),chipp(i,j),chipp(j,i), & + (alphasur(k,i,j),k=1,4),sigmap1(i,j),sigmap2(i,j), & + chis(i,j),chis(j,i), & + nstate(i,j),(wstate(k,i,j),k=1,4), & + dhead(1,1,i,j),dhead(1,2,i,j),dhead(2,1,i,j),dhead(2,2,i,j),& + dtail(1,i,j),dtail(2,i,j), & + epshead(i,j),sig0head(i,j), & + rborn(i,j),rborn(j,i),(wqdip(k,i,j),k=1,2),wquad(i,j), & + alphapol(i,j),alphapol(j,i), & + (alphiso(k,i,j),k=1,4),sigiso1(i,j),sigiso2(i,j),epsintab(i,j),debaykap(i,j) +! print *,eps(i,j),sigma(i,j),"SIGMAP",i,j,sigmap1(i,j),sigmap2(j,i) + END DO + END DO + DO i = 1, ntyp + DO j = i+1, ntyp + eps(i,j) = eps(j,i) + sigma(i,j) = sigma(j,i) + sigmap1(i,j)=sigmap1(j,i) + sigmap2(i,j)=sigmap2(j,i) + sigiso1(i,j)=sigiso1(j,i) + sigiso2(i,j)=sigiso2(j,i) +! print *,"ATU",sigma(j,i),sigma(i,j),i,j + nstate(i,j) = nstate(j,i) + dtail(1,i,j) = dtail(1,j,i) + dtail(2,i,j) = dtail(2,j,i) + DO k = 1, 4 + alphasur(k,i,j) = alphasur(k,j,i) + wstate(k,i,j) = wstate(k,j,i) + alphiso(k,i,j) = alphiso(k,j,i) + END DO + + dhead(2,1,i,j) = dhead(1,1,j,i) + dhead(2,2,i,j) = dhead(1,2,j,i) + dhead(1,1,i,j) = dhead(2,1,j,i) + dhead(1,2,i,j) = dhead(2,2,j,i) + + epshead(i,j) = epshead(j,i) + sig0head(i,j) = sig0head(j,i) + + DO k = 1, 2 + wqdip(k,i,j) = wqdip(k,j,i) + END DO + + wquad(i,j) = wquad(j,i) + epsintab(i,j) = epsintab(j,i) + debaykap(i,j)=debaykap(j,i) +! if (epsintab(i,j).ne.1.0) print *,"WHAT?",i,j,epsintab(i,j) + END DO + END DO + endif + ! goto 50 !--------------------- GBV potential ----------------------------------- case (5) @@ -1688,18 +2112,22 @@ enddo !el from module energy - COMMON.INTERACT------- ! allocate(aa(ntyp1,ntyp1),bb(ntyp1,ntyp1),chi(ntyp1,ntyp1)) !(ntyp,ntyp) - allocate(aa_aq(ntyp1,ntyp1),bb_aq(ntyp1,ntyp1),chi(ntyp1,ntyp1)) !(ntyp,ntyp) + if (.not.allocated(chi)) allocate(chi(ntyp1,ntyp1)) + allocate(aa_aq(ntyp1,ntyp1),bb_aq(ntyp1,ntyp1)) !(ntyp,ntyp) allocate(aa_lip(ntyp1,ntyp1),bb_lip(ntyp1,ntyp1)) !(ntyp,ntyp) - allocate(sigma(0:ntyp1,0:ntyp1),r0(ntyp1,ntyp1)) !(0:ntyp1,0:ntyp1) + if (.not.allocated(sigma)) allocate(sigma(0:ntyp1,0:ntyp1)) + allocate(r0(ntyp1,ntyp1)) !(0:ntyp1,0:ntyp1) do i=1,ntyp1 do j=1,ntyp1 aa_aq(i,j)=0.0D0 bb_aq(i,j)=0.0D0 aa_lip(i,j)=0.0d0 bb_lip(i,j)=0.0d0 + if (scelemode.eq.0) then chi(i,j)=0.0D0 sigma(i,j)=0.0D0 r0(i,j)=0.0D0 + endif enddo enddo !-------------------------------- @@ -1710,6 +2138,7 @@ enddo epslip(i,j)=epslip(j,i) enddo enddo + if (scelemode.eq.0) then do i=1,ntyp do j=i,ntyp sigma(i,j)=dsqrt(sigma0(i)**2+sigma0(j)**2) @@ -1718,6 +2147,7 @@ enddo rs0(j,i)=rs0(i,j) enddo enddo + endif if (lprint) write (iout,'(/a/10x,7a/72(1h-))') & 'Working parameters of the SC interactions:',& ' a ',' b ',' augm ',' sigma ',' r0 ',& @@ -1727,6 +2157,7 @@ enddo epsij=eps(i,j) if (ipot.eq.1 .or. ipot.eq.3 .or. ipot.eq.4) then rrij=sigma(i,j) + print *,"rrij",rrij else rrij=rr0(i)+rr0(j) endif @@ -1747,7 +2178,7 @@ enddo bb_lip(i,j)=-sigeps*epsijlip*rrij aa_lip(j,i)=aa_lip(i,j) bb_lip(j,i)=bb_lip(i,j) - if (ipot.gt.2) then + if ((ipot.gt.2).and. (scelemode.eq.0))then sigt1sq=sigma0(i)**2 sigt2sq=sigma0(j)**2 sigii1=sigii(i) @@ -1952,9 +2383,12 @@ enddo !----------------------------------------------------------------------------- subroutine read_general_data(*) - use control_data, only:indpdb,symetr,r_cut_ele,rlamb_ele,ions - use energy_data, only:distchainmax - use geometry_data, only:boxxsize,boxysize,boxzsize + use control_data, only:indpdb,symetr,r_cut_ele,rlamb_ele,ions,& + scelemode,TUBEmode,tor_mode + + use energy_data, only:distchainmax,tubeR0,tubecenter + use geometry_data, only:boxxsize,boxysize,boxzsize,bordtubetop,& + bordtubebot,tubebufthick,buftubebot,buftubetop ! implicit none ! include "DIMENSIONS" ! include "DIMENSIONS.ZSCOPT" @@ -2034,6 +2468,25 @@ enddo call reada(controlcard,'BOXX',boxxsize,100.0d0) call reada(controlcard,'BOXY',boxysize,100.0d0) call reada(controlcard,'BOXZ',boxzsize,100.0d0) + call readi(controlcard,"SCELEMODE",scelemode,0) + print *,"SCELE",scelemode + call readi(controlcard,'TORMODE',tor_mode,0) +!C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then + write(iout,*) "torsional and valence angle mode",tor_mode + + call readi(controlcard,'TUBEMOD',tubemode,0) + + if (TUBEmode.gt.0) then + call reada(controlcard,"XTUBE",tubecenter(1),0.0d0) + call reada(controlcard,"YTUBE",tubecenter(2),0.0d0) + call reada(controlcard,"ZTUBE",tubecenter(3),0.0d0) + call reada(controlcard,"RTUBE",tubeR0,0.0d0) + call reada(controlcard,"TUBETOP",bordtubetop,boxzsize) + call reada(controlcard,"TUBEBOT",bordtubebot,0.0d0) + call reada(controlcard,"TUBEBUF",tubebufthick,1.0d0) + buftubebot=bordtubebot+tubebufthick + buftubetop=bordtubetop-tubebufthick + endif ions=index(controlcard,"IONS").gt.0 call reada(controlcard,"R_CUT_ELE",r_cut_ele,15.0d0) call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0) diff --git a/source/wham/wham_calc.F90 b/source/wham/wham_calc.F90 index 4182897..d1ed2f7 100644 --- a/source/wham/wham_calc.F90 +++ b/source/wham/wham_calc.F90 @@ -214,7 +214,7 @@ MPI_MIN,WHAM_COMM,IERROR) call MPI_AllReduce(rgymax,rgymax_t,1,MPI_DOUBLE_PRECISION,& MPI_MAX,WHAM_COMM,IERROR) - potEmin=potEmin_t/2 + potEmin=potEmin_t !/2 try now?? rgymin=rgymin_t rgymax=rgymax_t rmsmin=rmsmin_t @@ -253,7 +253,7 @@ !#ifdef DEBUG write (iout,*) wsc,wscp,welec,wvdwpp,wang,wtor,wscloc,& wcorr,wcorr5,wcorr6,wturn4,wturn3,wturn6,wel_loc,& - wtor_d,wsccor,wbond + wtor_d,wsccor,wbond,wcatcat !#endif do ib=1,nT_h(iparm) !el old rascale weights @@ -343,12 +343,12 @@ esccor=enetb(21,i,iparm) ! edihcnstr=enetb(20,i,iparm) edihcnstr=enetb(19,i,iparm) - ecationcation=enetb(42,i,iparm) - ecation_prot=enetb(41,i,iparm) + ecationcation=enetb(41,i,iparm) + ecation_prot=enetb(42,i,iparm) #ifdef DEBUG - write (iout,'(3i5,6f5.2,14f12.3)') i,ib,iparm,(ft(l),l=1,6),& + write (iout,'(3i5,6f5.2,15f12.3)') i,ib,iparm,(ft(l),l=1,6),& evdw+evdw_t,evdw2,ees,evdw1,ecorr,eel_loc,estr,ebe,escloc,& - etors,etors_d,eello_turn3,eello_turn4,esccor + etors,etors_d,eello_turn3,eello_turn4,esccor,ecationcation #endif !#ifdef SPLITELE @@ -755,8 +755,8 @@ ! edihcnstr=enetb(20,t,iparm) edihcnstr=enetb(19,t,iparm) edihcnstr=0.0d0 - ecationcation=enetb(42,t,iparm) - ecation_prot=enetb(41,t,iparm) + ecationcation=enetb(41,t,iparm) + ecation_prot=enetb(42,t,iparm) do k=0,nGridT betaT=startGridT+k*delta_T -- 1.7.9.5