diff --git a/CHANGES/v2.5.md b/CHANGES/v2.5.md index 840af4453a69a44f8f41b4e8070c9c3def39d51b..43541fb60104e7aeb03cb73d6e3bccbf0874932e 100644 --- a/CHANGES/v2.5.md +++ b/CHANGES/v2.5.md @@ -44,6 +44,12 @@ Changes from version 2.4 which are relevant for users: - \ref completion (used to generate command line completion scripts). - \ref pdbrenumber (see \issue{371}). +- Changes in the ISDB module + - \ref CS2BACKBONE is now mpi parallelised in particular with DOSCORE and CAMSHIFT + - \ref SAXS has an additional implementation based on bessel functions that can be faster for large systems (new keyword BESSEL) + - \ref METAINFERENCE and all related methods has a new keyword REGRES_ZERO to scale data using a linear scale fit + - \ref CALIBER new bias to perform Maximum Caliber replica-averaged restrained simulations + - Changes in the eABF/DRR module (contributed by Haochuan Chen and Haohao Fu): - \ref DRR now supports the extended generalized ABF(egABF) method. - \ref DRR accepts different GRID options for CVs and extended variables. @@ -59,8 +65,6 @@ Changes from version 2.4 which are relevant for users: - Other changes: - \ref EXTERNAL can now SCALE the input grid. This allows for more flexibility without modifying the grid file. - \ref ALPHABETA can now combine dihedrals with different coefficients - - \ref CS2BACKBONE (isdb module) is now mpi parallelised in particular with DOSCORE and CAMSHIFT - - \ref SAXS (isdb module) there is an additional implementation based on bessel functions that can be faster for large systems (new keyword BESSEL) - \ref INCLUDE can now be used also before setup actions. - Libmatheval is not used anymore. \ref MATHEVAL (and \ref CUSTOM) are still available but employ an internal implementation of the lepton library. diff --git a/regtest/isdb/rt-caliber/Makefile b/regtest/isdb/rt-caliber/Makefile new file mode 100644 index 0000000000000000000000000000000000000000..3703b27cea227aa053fb6d1d73f861e4384dbcee --- /dev/null +++ b/regtest/isdb/rt-caliber/Makefile @@ -0,0 +1 @@ +include ../../scripts/test.make diff --git a/regtest/isdb/rt-caliber/Qaverages.dat b/regtest/isdb/rt-caliber/Qaverages.dat new file mode 100644 index 0000000000000000000000000000000000000000..f37ddfc9eff1bfaa03eff6ad87f62bf9f0b2dc33 --- /dev/null +++ b/regtest/isdb/rt-caliber/Qaverages.dat @@ -0,0 +1,501 @@ +0.000000000000000000e+00 7.443500000000052852e-02 +1.000000000000000000e+03 4.782885399999999010e-02 +2.000000000000000000e+03 5.607749799999997564e-02 +3.000000000000000000e+03 6.611566999999998751e-02 +4.000000000000000000e+03 6.298131599999996777e-02 +5.000000000000000000e+03 5.650602000000000402e-02 +6.000000000000000000e+03 5.914363799999997756e-02 +7.000000000000000000e+03 6.399480400000003011e-02 +8.000000000000000000e+03 6.554257400000007561e-02 +9.000000000000000000e+03 7.244970400000001753e-02 +1.000000000000000000e+04 8.297708200000010470e-02 +1.100000000000000000e+04 9.649917999999993445e-02 +1.200000000000000000e+04 1.027288260000000369e-01 +1.300000000000000000e+04 1.057008019999999965e-01 +1.400000000000000000e+04 1.079640280000000035e-01 +1.500000000000000000e+04 1.117816940000000703e-01 +1.600000000000000000e+04 1.108070040000001255e-01 +1.700000000000000000e+04 1.134541899999999687e-01 +1.800000000000000000e+04 1.098611400000000515e-01 +1.900000000000000000e+04 1.072102019999999906e-01 +2.000000000000000000e+04 1.085445980000000615e-01 +2.100000000000000000e+04 1.135731060000000348e-01 +2.200000000000000000e+04 1.191038719999999024e-01 +2.300000000000000000e+04 1.243017319999999565e-01 +2.400000000000000000e+04 1.201812720000000750e-01 +2.500000000000000000e+04 1.157263260000000599e-01 +2.600000000000000000e+04 1.100925539999999814e-01 +2.700000000000000000e+04 1.109309640000000485e-01 +2.800000000000000000e+04 1.121387340000000871e-01 +2.900000000000000000e+04 1.158201319999999784e-01 +3.000000000000000000e+04 1.178793520000000067e-01 +3.100000000000000000e+04 1.190168120000000135e-01 +3.200000000000000000e+04 1.205617380000000016e-01 +3.300000000000000000e+04 1.243166840000000384e-01 +3.400000000000000000e+04 1.266460619999998649e-01 +3.500000000000000000e+04 1.320553279999999441e-01 +3.600000000000000000e+04 1.343793439999999284e-01 +3.700000000000000000e+04 1.299888600000000116e-01 +3.800000000000000000e+04 1.281730420000000703e-01 +3.900000000000000000e+04 1.312398819999999744e-01 +4.000000000000000000e+04 1.348972639999999890e-01 +4.100000000000000000e+04 1.390801140000000602e-01 +4.200000000000000000e+04 1.402975500000001874e-01 +4.300000000000000000e+04 1.384006460000000716e-01 +4.400000000000000000e+04 1.358499020000000224e-01 +4.500000000000000000e+04 1.356443400000000021e-01 +4.600000000000000000e+04 1.367546680000000514e-01 +4.700000000000000000e+04 1.428825939999999739e-01 +4.800000000000000000e+04 1.437093060000000089e-01 +4.900000000000000000e+04 1.411044660000000950e-01 +5.000000000000000000e+04 1.393018699999998555e-01 +5.100000000000000000e+04 1.392468059999998897e-01 +5.200000000000000000e+04 1.427314780000001615e-01 +5.300000000000000000e+04 1.494656539999998646e-01 +5.400000000000000000e+04 1.544176759999999760e-01 +5.500000000000000000e+04 1.601531980000000521e-01 +5.600000000000000000e+04 1.642821020000000687e-01 +5.700000000000000000e+04 1.665094280000001814e-01 +5.800000000000000000e+04 1.656651459999998466e-01 +5.900000000000000000e+04 1.648939060000001178e-01 +6.000000000000000000e+04 1.612244780000000322e-01 +6.100000000000000000e+04 1.593273460000000918e-01 +6.200000000000000000e+04 1.592322460000000495e-01 +6.300000000000000000e+04 1.618426379999999831e-01 +6.400000000000000000e+04 1.630193560000001707e-01 +6.500000000000000000e+04 1.678868259999998502e-01 +6.600000000000000000e+04 1.712108079999998533e-01 +6.700000000000000000e+04 1.826008439999999289e-01 +6.800000000000000000e+04 1.909277299999998512e-01 +6.900000000000000000e+04 1.952022999999999120e-01 +7.000000000000000000e+04 2.022678760000000131e-01 +7.100000000000000000e+04 2.061609599999998654e-01 +7.200000000000000000e+04 2.082702600000000959e-01 +7.300000000000000000e+04 2.070076119999999797e-01 +7.400000000000000000e+04 2.040152020000001731e-01 +7.500000000000000000e+04 2.015203700000000875e-01 +7.600000000000000000e+04 2.005693100000000006e-01 +7.700000000000000000e+04 1.973828099999997976e-01 +7.800000000000000000e+04 1.920514039999999256e-01 +7.900000000000000000e+04 1.932291359999999958e-01 +8.000000000000000000e+04 1.928419620000001999e-01 +8.100000000000000000e+04 1.919228440000000924e-01 +8.200000000000000000e+04 1.879411080000000234e-01 +8.300000000000000000e+04 1.889222140000001993e-01 +8.400000000000000000e+04 1.908450799999998893e-01 +8.500000000000000000e+04 1.879246360000000615e-01 +8.600000000000000000e+04 1.833099420000000035e-01 +8.700000000000000000e+04 1.817689580000000915e-01 +8.800000000000000000e+04 1.809812700000000829e-01 +8.900000000000000000e+04 1.800451880000000504e-01 +9.000000000000000000e+04 1.811042999999999126e-01 +9.100000000000000000e+04 1.819289280000001008e-01 +9.200000000000000000e+04 1.846117339999997775e-01 +9.300000000000000000e+04 1.894187079999999912e-01 +9.400000000000000000e+04 1.851843340000000337e-01 +9.500000000000000000e+04 1.851542020000000177e-01 +9.600000000000000000e+04 1.847731340000000333e-01 +9.700000000000000000e+04 1.856753539999999592e-01 +9.800000000000000000e+04 1.894862740000000934e-01 +9.900000000000000000e+04 1.938538880000001130e-01 +1.000000000000000000e+05 2.002189120000000266e-01 +1.010000000000000000e+05 2.079729459999997143e-01 +1.020000000000000000e+05 2.151246519999999995e-01 +1.030000000000000000e+05 2.222804540000000439e-01 +1.040000000000000000e+05 2.282939099999998223e-01 +1.050000000000000000e+05 2.346418879999999652e-01 +1.060000000000000000e+05 2.340413640000001960e-01 +1.070000000000000000e+05 2.387018720000000371e-01 +1.080000000000000000e+05 2.446861220000001724e-01 +1.090000000000000000e+05 2.506001459999998016e-01 +1.100000000000000000e+05 2.584995220000002036e-01 +1.110000000000000000e+05 2.639496380000000419e-01 +1.120000000000000000e+05 2.660931860000000371e-01 +1.130000000000000000e+05 2.694745159999999973e-01 +1.140000000000000000e+05 2.713951019999999992e-01 +1.150000000000000000e+05 2.711399359999998038e-01 +1.160000000000000000e+05 2.729090380000003258e-01 +1.170000000000000000e+05 2.742953220000002856e-01 +1.180000000000000000e+05 2.745139280000000181e-01 +1.190000000000000000e+05 2.758631840000000945e-01 +1.200000000000000000e+05 2.784365099999999149e-01 +1.210000000000000000e+05 2.788716400000000872e-01 +1.220000000000000000e+05 2.795883699999998084e-01 +1.230000000000000000e+05 2.812095340000002608e-01 +1.240000000000000000e+05 2.834075520000000070e-01 +1.250000000000000000e+05 2.865440560000000469e-01 +1.260000000000000000e+05 2.874799340000001036e-01 +1.270000000000000000e+05 2.916171800000000869e-01 +1.280000000000000000e+05 2.968801619999998920e-01 +1.290000000000000000e+05 3.035590879999996994e-01 +1.300000000000000000e+05 3.100981859999999424e-01 +1.310000000000000000e+05 3.181280160000000135e-01 +1.320000000000000000e+05 3.231979880000000471e-01 +1.330000000000000000e+05 3.278055980000000869e-01 +1.340000000000000000e+05 3.320442620000005896e-01 +1.350000000000000000e+05 3.345744219999998714e-01 +1.360000000000000000e+05 3.370458259999998263e-01 +1.370000000000000000e+05 3.395394619999998476e-01 +1.380000000000000000e+05 3.422543899999998529e-01 +1.390000000000000000e+05 3.477584679999998762e-01 +1.400000000000000000e+05 3.488078419999999236e-01 +1.410000000000000000e+05 3.501677639999997704e-01 +1.420000000000000000e+05 3.541321000000001162e-01 +1.430000000000000000e+05 3.581256880000001086e-01 +1.440000000000000000e+05 3.633888520000003175e-01 +1.450000000000000000e+05 3.668034299999995973e-01 +1.460000000000000000e+05 3.683116280000003351e-01 +1.470000000000000000e+05 3.728373040000002581e-01 +1.480000000000000000e+05 3.750480519999999096e-01 +1.490000000000000000e+05 3.789673120000000561e-01 +1.500000000000000000e+05 3.851196719999995244e-01 +1.510000000000000000e+05 3.910460439999998150e-01 +1.520000000000000000e+05 3.928223919999997427e-01 +1.530000000000000000e+05 3.953805139999997387e-01 +1.540000000000000000e+05 3.967094200000001458e-01 +1.550000000000000000e+05 3.980394360000000242e-01 +1.560000000000000000e+05 4.014043779999999506e-01 +1.570000000000000000e+05 4.039606019999995579e-01 +1.580000000000000000e+05 4.067817699999997094e-01 +1.590000000000000000e+05 4.105731600000002701e-01 +1.600000000000000000e+05 4.155473119999997800e-01 +1.610000000000000000e+05 4.174706539999998300e-01 +1.620000000000000000e+05 4.228700780000000381e-01 +1.630000000000000000e+05 4.251336719999999625e-01 +1.640000000000000000e+05 4.280854760000001313e-01 +1.650000000000000000e+05 4.337042080000000910e-01 +1.660000000000000000e+05 4.406500480000002384e-01 +1.670000000000000000e+05 4.433785819999998545e-01 +1.680000000000000000e+05 4.482019460000001287e-01 +1.690000000000000000e+05 4.522862300000006219e-01 +1.700000000000000000e+05 4.550641279999997346e-01 +1.710000000000000000e+05 4.559146360000003173e-01 +1.720000000000000000e+05 4.592384980000001060e-01 +1.730000000000000000e+05 4.631929039999999609e-01 +1.740000000000000000e+05 4.666824580000000222e-01 +1.750000000000000000e+05 4.715138359999998530e-01 +1.760000000000000000e+05 4.794720300000000490e-01 +1.770000000000000000e+05 4.834903360000000205e-01 +1.780000000000000000e+05 4.844139159999997779e-01 +1.790000000000000000e+05 4.846910160000003631e-01 +1.800000000000000000e+05 4.847626020000000979e-01 +1.810000000000000000e+05 4.854316560000003222e-01 +1.820000000000000000e+05 4.893824339999998663e-01 +1.830000000000000000e+05 4.906448800000004495e-01 +1.840000000000000000e+05 4.922080160000001925e-01 +1.850000000000000000e+05 4.946878179999998348e-01 +1.860000000000000000e+05 4.987338980000000088e-01 +1.870000000000000000e+05 5.022803780000005824e-01 +1.880000000000000000e+05 5.034389820000000348e-01 +1.890000000000000000e+05 5.031699379999997612e-01 +1.900000000000000000e+05 5.049028639999998402e-01 +1.910000000000000000e+05 5.099838360000003012e-01 +1.920000000000000000e+05 5.157410580000006961e-01 +1.930000000000000000e+05 5.182363039999997589e-01 +1.940000000000000000e+05 5.197325480000002385e-01 +1.950000000000000000e+05 5.223613799999999863e-01 +1.960000000000000000e+05 5.256032359999999182e-01 +1.970000000000000000e+05 5.314073419999998382e-01 +1.980000000000000000e+05 5.373958780000003266e-01 +1.990000000000000000e+05 5.440329399999997984e-01 +2.000000000000000000e+05 5.488316080000009700e-01 +2.010000000000000000e+05 5.503979199999999850e-01 +2.020000000000000000e+05 5.529269720000002941e-01 +2.030000000000000000e+05 5.571756460000001354e-01 +2.040000000000000000e+05 5.636313860000001785e-01 +2.050000000000000000e+05 5.679959860000005634e-01 +2.060000000000000000e+05 5.708596280000003409e-01 +2.070000000000000000e+05 5.741887179999994872e-01 +2.080000000000000000e+05 5.780932680000002710e-01 +2.090000000000000000e+05 5.826644119999998539e-01 +2.100000000000000000e+05 5.861513480000003495e-01 +2.110000000000000000e+05 5.902139320000008293e-01 +2.120000000000000000e+05 5.938952760000003606e-01 +2.130000000000000000e+05 5.987347040000003950e-01 +2.140000000000000000e+05 6.032528600000003349e-01 +2.150000000000000000e+05 6.078451939999994780e-01 +2.160000000000000000e+05 6.095975019999998468e-01 +2.170000000000000000e+05 6.136406220000009126e-01 +2.180000000000000000e+05 6.156279300000002674e-01 +2.190000000000000000e+05 6.178962520000002012e-01 +2.200000000000000000e+05 6.206229240000001868e-01 +2.210000000000000000e+05 6.228219680000005587e-01 +2.220000000000000000e+05 6.255704740000006536e-01 +2.230000000000000000e+05 6.279148600000000746e-01 +2.240000000000000000e+05 6.310831980000005670e-01 +2.250000000000000000e+05 6.345840660000003908e-01 +2.260000000000000000e+05 6.392222240000010602e-01 +2.270000000000000000e+05 6.436906259999995994e-01 +2.280000000000000000e+05 6.484052120000006747e-01 +2.290000000000000000e+05 6.536483540000003201e-01 +2.300000000000000000e+05 6.578628939999998648e-01 +2.310000000000000000e+05 6.627671679999999066e-01 +2.320000000000000000e+05 6.662007220000004670e-01 +2.330000000000000000e+05 6.700120600000003535e-01 +2.340000000000000000e+05 6.736105280000000972e-01 +2.350000000000000000e+05 6.775587440000004902e-01 +2.360000000000000000e+05 6.804812579999999222e-01 +2.370000000000000000e+05 6.822766159999996693e-01 +2.380000000000000000e+05 6.838889860000005871e-01 +2.390000000000000000e+05 6.871059139999996512e-01 +2.400000000000000000e+05 6.906701040000003955e-01 +2.410000000000000000e+05 6.969089119999997139e-01 +2.420000000000000000e+05 7.014348999999998613e-01 +2.430000000000000000e+05 7.062068820000007019e-01 +2.440000000000000000e+05 7.086591540000001244e-01 +2.450000000000000000e+05 7.096774499999999319e-01 +2.460000000000000000e+05 7.118449840000002915e-01 +2.470000000000000000e+05 7.136841760000004475e-01 +2.480000000000000000e+05 7.134710540000007706e-01 +2.490000000000000000e+05 7.150336640000000132e-01 +2.500000000000000000e+05 7.176564380000010068e-01 +2.510000000000000000e+05 7.202942220000003726e-01 +2.520000000000000000e+05 7.248647679999995486e-01 +2.530000000000000000e+05 7.293641420000003261e-01 +2.540000000000000000e+05 7.323074619999999646e-01 +2.550000000000000000e+05 7.347425180000002332e-01 +2.560000000000000000e+05 7.374476820000002153e-01 +2.570000000000000000e+05 7.402872100000000843e-01 +2.580000000000000000e+05 7.439755560000007595e-01 +2.590000000000000000e+05 7.469396300000007427e-01 +2.600000000000000000e+05 7.498830140000009026e-01 +2.610000000000000000e+05 7.542468200000007617e-01 +2.620000000000000000e+05 7.582247940000007569e-01 +2.630000000000000000e+05 7.626353400000011051e-01 +2.640000000000000000e+05 7.671383420000003062e-01 +2.650000000000000000e+05 7.690333860000008182e-01 +2.660000000000000000e+05 7.727548980000008561e-01 +2.670000000000000000e+05 7.750374200000006430e-01 +2.680000000000000000e+05 7.760245440000006489e-01 +2.690000000000000000e+05 7.780812879999999820e-01 +2.700000000000000000e+05 7.831540080000002346e-01 +2.710000000000000000e+05 7.870012700000003081e-01 +2.720000000000000000e+05 7.896707500000008650e-01 +2.730000000000000000e+05 7.923290180000002181e-01 +2.740000000000000000e+05 7.972646560000000715e-01 +2.750000000000000000e+05 7.976267580000007129e-01 +2.760000000000000000e+05 7.996154880000007070e-01 +2.770000000000000000e+05 8.023490760000008271e-01 +2.780000000000000000e+05 8.052161560000008489e-01 +2.790000000000000000e+05 8.064512100000006400e-01 +2.800000000000000000e+05 8.090262020000005272e-01 +2.810000000000000000e+05 8.125600680000002463e-01 +2.820000000000000000e+05 8.142756900000005515e-01 +2.830000000000000000e+05 8.160210260000008153e-01 +2.840000000000000000e+05 8.177618080000004230e-01 +2.850000000000000000e+05 8.212098860000006662e-01 +2.860000000000000000e+05 8.235971380000005615e-01 +2.870000000000000000e+05 8.268775120000013690e-01 +2.880000000000000000e+05 8.290138920000006406e-01 +2.890000000000000000e+05 8.315856600000011145e-01 +2.900000000000000000e+05 8.331459540000004660e-01 +2.910000000000000000e+05 8.337720959999999071e-01 +2.920000000000000000e+05 8.355878520000006304e-01 +2.930000000000000000e+05 8.377243580000005574e-01 +2.940000000000000000e+05 8.391014760000004014e-01 +2.950000000000000000e+05 8.414871680000003407e-01 +2.960000000000000000e+05 8.427632819999995029e-01 +2.970000000000000000e+05 8.448468120000004467e-01 +2.980000000000000000e+05 8.467773460000000973e-01 +2.990000000000000000e+05 8.496589220000000386e-01 +3.000000000000000000e+05 8.505874340000006972e-01 +3.010000000000000000e+05 8.518655440000005852e-01 +3.020000000000000000e+05 8.542549639999997835e-01 +3.030000000000000000e+05 8.546543140000008032e-01 +3.040000000000000000e+05 8.555541620000011171e-01 +3.050000000000000000e+05 8.573298180000012980e-01 +3.060000000000000000e+05 8.606145240000002694e-01 +3.070000000000000000e+05 8.618253120000008982e-01 +3.080000000000000000e+05 8.628441340000011506e-01 +3.090000000000000000e+05 8.631881660000015888e-01 +3.100000000000000000e+05 8.646507860000007817e-01 +3.110000000000000000e+05 8.660334420000008748e-01 +3.120000000000000000e+05 8.670479200000006381e-01 +3.130000000000000000e+05 8.688061960000000861e-01 +3.140000000000000000e+05 8.712195060000000879e-01 +3.150000000000000000e+05 8.718409180000011860e-01 +3.160000000000000000e+05 8.715725580000012185e-01 +3.170000000000000000e+05 8.720319139999995750e-01 +3.180000000000000000e+05 8.737934520000004701e-01 +3.190000000000000000e+05 8.747422159999993507e-01 +3.200000000000000000e+05 8.751917680000007582e-01 +3.210000000000000000e+05 8.756721739999993037e-01 +3.220000000000000000e+05 8.765817360000003333e-01 +3.230000000000000000e+05 8.775885900000000017e-01 +3.240000000000000000e+05 8.794875080000004175e-01 +3.250000000000000000e+05 8.802486900000001668e-01 +3.260000000000000000e+05 8.821218420000014326e-01 +3.270000000000000000e+05 8.829647660000017328e-01 +3.280000000000000000e+05 8.837307800000012703e-01 +3.290000000000000000e+05 8.851845940000023782e-01 +3.300000000000000000e+05 8.874981000000008446e-01 +3.310000000000000000e+05 8.869257319999999112e-01 +3.320000000000000000e+05 8.869456180000002687e-01 +3.330000000000000000e+05 8.881574359999999668e-01 +3.340000000000000000e+05 8.897856680000010288e-01 +3.350000000000000000e+05 8.914885220000006161e-01 +3.360000000000000000e+05 8.929916180000012638e-01 +3.370000000000000000e+05 8.948782960000007947e-01 +3.380000000000000000e+05 8.943573840000005326e-01 +3.390000000000000000e+05 8.953558260000008540e-01 +3.400000000000000000e+05 8.966014620000004598e-01 +3.410000000000000000e+05 8.984813659999997260e-01 +3.420000000000000000e+05 8.983560220000005314e-01 +3.430000000000000000e+05 8.987468600000007024e-01 +3.440000000000000000e+05 8.993468800000009589e-01 +3.450000000000000000e+05 8.998362540000005572e-01 +3.460000000000000000e+05 9.007040100000003324e-01 +3.470000000000000000e+05 9.016307020000008388e-01 +3.480000000000000000e+05 9.022230819999998985e-01 +3.490000000000000000e+05 9.016510500000005646e-01 +3.500000000000000000e+05 9.010675400000007773e-01 +3.510000000000000000e+05 9.005514340000013718e-01 +3.520000000000000000e+05 9.012065579999998244e-01 +3.530000000000000000e+05 9.036295739999997689e-01 +3.540000000000000000e+05 9.051161460000007875e-01 +3.550000000000000000e+05 9.065767359999999941e-01 +3.560000000000000000e+05 9.080495680000004732e-01 +3.570000000000000000e+05 9.093111840000008828e-01 +3.580000000000000000e+05 9.108024720000008623e-01 +3.590000000000000000e+05 9.109734980000012694e-01 +3.600000000000000000e+05 9.103190000000010995e-01 +3.610000000000000000e+05 9.115989820000006683e-01 +3.620000000000000000e+05 9.125968380000006874e-01 +3.630000000000000000e+05 9.131023660000007203e-01 +3.640000000000000000e+05 9.148243260000006316e-01 +3.650000000000000000e+05 9.139154899999996085e-01 +3.660000000000000000e+05 9.133511940000011986e-01 +3.670000000000000000e+05 9.137394820000006579e-01 +3.680000000000000000e+05 9.148641140000000060e-01 +3.690000000000000000e+05 9.160044900000005041e-01 +3.700000000000000000e+05 9.166018300000006169e-01 +3.710000000000000000e+05 9.171233760000011559e-01 +3.720000000000000000e+05 9.178556779999996751e-01 +3.730000000000000000e+05 9.178648420000002917e-01 +3.740000000000000000e+05 9.187414300000013023e-01 +3.750000000000000000e+05 9.205547339999993195e-01 +3.760000000000000000e+05 9.220803260000012269e-01 +3.770000000000000000e+05 9.206663020000002140e-01 +3.780000000000000000e+05 9.196974240000006517e-01 +3.790000000000000000e+05 9.214880280000012913e-01 +3.800000000000000000e+05 9.225715640000012874e-01 +3.810000000000000000e+05 9.233609720000004462e-01 +3.820000000000000000e+05 9.246728340000002211e-01 +3.830000000000000000e+05 9.245971480000008125e-01 +3.840000000000000000e+05 9.245765820000000357e-01 +3.850000000000000000e+05 9.241937800000007419e-01 +3.860000000000000000e+05 9.243753080000004507e-01 +3.870000000000000000e+05 9.246531920000004012e-01 +3.880000000000000000e+05 9.263928519999999889e-01 +3.890000000000000000e+05 9.267628020000009403e-01 +3.900000000000000000e+05 9.261656220000010764e-01 +3.910000000000000000e+05 9.268191400000010960e-01 +3.920000000000000000e+05 9.266219780000003459e-01 +3.930000000000000000e+05 9.280437040000000248e-01 +3.940000000000000000e+05 9.285165260000007859e-01 +3.950000000000000000e+05 9.283971500000003374e-01 +3.960000000000000000e+05 9.284482420000006453e-01 +3.970000000000000000e+05 9.294999540000011828e-01 +3.980000000000000000e+05 9.304906760000012378e-01 +3.990000000000000000e+05 9.302026940000006627e-01 +4.000000000000000000e+05 9.315608640000004881e-01 +4.010000000000000000e+05 9.321296280000022083e-01 +4.020000000000000000e+05 9.330253600000008030e-01 +4.030000000000000000e+05 9.333381220000007694e-01 +4.040000000000000000e+05 9.340813240000007678e-01 +4.050000000000000000e+05 9.340758540000007049e-01 +4.060000000000000000e+05 9.334881020000014296e-01 +4.070000000000000000e+05 9.348207500000005332e-01 +4.080000000000000000e+05 9.355898800000007620e-01 +4.090000000000000000e+05 9.368991819999999970e-01 +4.100000000000000000e+05 9.364489820000008180e-01 +4.110000000000000000e+05 9.370557240000009225e-01 +4.120000000000000000e+05 9.375718480000003741e-01 +4.130000000000000000e+05 9.373707680000015197e-01 +4.140000000000000000e+05 9.371115740000003358e-01 +4.150000000000000000e+05 9.373151100000014502e-01 +4.160000000000000000e+05 9.367829600000010526e-01 +4.170000000000000000e+05 9.368473660000010694e-01 +4.180000000000000000e+05 9.365574760000012766e-01 +4.190000000000000000e+05 9.374192500000003703e-01 +4.200000000000000000e+05 9.389462440000012622e-01 +4.210000000000000000e+05 9.401917279999996158e-01 +4.220000000000000000e+05 9.422190740000006004e-01 +4.230000000000000000e+05 9.432876280000007663e-01 +4.240000000000000000e+05 9.446560100000003235e-01 +4.250000000000000000e+05 9.453807600000001532e-01 +4.260000000000000000e+05 9.446751180000011461e-01 +4.270000000000000000e+05 9.457434340000011597e-01 +4.280000000000000000e+05 9.465151360000003677e-01 +4.290000000000000000e+05 9.472960860000011474e-01 +4.300000000000000000e+05 9.482765000000010769e-01 +4.310000000000000000e+05 9.478908800000005463e-01 +4.320000000000000000e+05 9.460821020000003134e-01 +4.330000000000000000e+05 9.448400760000009724e-01 +4.340000000000000000e+05 9.445204700000006115e-01 +4.350000000000000000e+05 9.461976080000005229e-01 +4.360000000000000000e+05 9.475175080000012295e-01 +4.370000000000000000e+05 9.474042760000006558e-01 +4.380000000000000000e+05 9.454706860000006152e-01 +4.390000000000000000e+05 9.444982940000001825e-01 +4.400000000000000000e+05 9.448987039999997561e-01 +4.410000000000000000e+05 9.457109180000007331e-01 +4.420000000000000000e+05 9.467356160000014453e-01 +4.430000000000000000e+05 9.467826460000009492e-01 +4.440000000000000000e+05 9.479620000000010815e-01 +4.450000000000000000e+05 9.474801980000011614e-01 +4.460000000000000000e+05 9.481217580000005363e-01 +4.470000000000000000e+05 9.494399399999997602e-01 +4.480000000000000000e+05 9.495412160000004098e-01 +4.490000000000000000e+05 9.496306179999999264e-01 +4.500000000000000000e+05 9.494177180000010763e-01 +4.510000000000000000e+05 9.504949560000007791e-01 +4.520000000000000000e+05 9.503106500000008117e-01 +4.530000000000000000e+05 9.495979660000011258e-01 +4.540000000000000000e+05 9.484641540000015025e-01 +4.550000000000000000e+05 9.476713520000007440e-01 +4.560000000000000000e+05 9.485218100000020192e-01 +4.570000000000000000e+05 9.484525980000014789e-01 +4.580000000000000000e+05 9.484553900000017856e-01 +4.590000000000000000e+05 9.499923020000007323e-01 +4.600000000000000000e+05 9.509278379999998032e-01 +4.610000000000000000e+05 9.510344320000014839e-01 +4.620000000000000000e+05 9.510952000000003626e-01 +4.630000000000000000e+05 9.518688100000003427e-01 +4.640000000000000000e+05 9.532756260000010418e-01 +4.650000000000000000e+05 9.537254180000011017e-01 +4.660000000000000000e+05 9.531861600000006707e-01 +4.670000000000000000e+05 9.528920680000005916e-01 +4.680000000000000000e+05 9.514717220000005193e-01 +4.690000000000000000e+05 9.510398460000012433e-01 +4.700000000000000000e+05 9.496716620000004161e-01 +4.710000000000000000e+05 9.505725479999999061e-01 +4.720000000000000000e+05 9.510045580000013876e-01 +4.730000000000000000e+05 9.513371700000015085e-01 +4.740000000000000000e+05 9.520925060000011708e-01 +4.750000000000000000e+05 9.528947900000011018e-01 +4.760000000000000000e+05 9.527411940000009016e-01 +4.770000000000000000e+05 9.525942320000012620e-01 +4.780000000000000000e+05 9.529587940000008306e-01 +4.790000000000000000e+05 9.522571580000007696e-01 +4.800000000000000000e+05 9.519584440000010694e-01 +4.810000000000000000e+05 9.521757620000008693e-01 +4.820000000000000000e+05 9.517025020000009494e-01 +4.830000000000000000e+05 9.506951140000011735e-01 +4.840000000000000000e+05 9.506895140000003463e-01 +4.850000000000000000e+05 9.515789480000013123e-01 +4.860000000000000000e+05 9.529097580000004397e-01 +4.870000000000000000e+05 9.546223500000001749e-01 +4.880000000000000000e+05 9.538235020000011000e-01 +4.890000000000000000e+05 9.540642640000003283e-01 +4.900000000000000000e+05 9.543749980000004740e-01 +4.910000000000000000e+05 9.544230020000015191e-01 +4.920000000000000000e+05 9.533261540000006473e-01 +4.930000000000000000e+05 9.530456119999999309e-01 +4.940000000000000000e+05 9.532780400000012699e-01 +4.950000000000000000e+05 9.532065020000006772e-01 +4.960000000000000000e+05 9.535974080000005904e-01 +4.970000000000000000e+05 9.531710240000005330e-01 +4.980000000000000000e+05 9.534724620000011308e-01 +4.990000000000000000e+05 9.539503000000000288e-01 +5.000000000000000000e+05 9.549465360000011227e-01 diff --git a/regtest/isdb/rt-caliber/config b/regtest/isdb/rt-caliber/config new file mode 100644 index 0000000000000000000000000000000000000000..41685b007b88843b6f7bb81b26ea9b19bd9bd986 --- /dev/null +++ b/regtest/isdb/rt-caliber/config @@ -0,0 +1,3 @@ +mpiprocs=4 +type=driver +arg="--plumed plumed.dat --trajectory-stride 1 --timestep 0.005 --mf_xtc trajout.xtc --multi 4 --dump-forces allforce" diff --git a/regtest/isdb/rt-caliber/first.pdb b/regtest/isdb/rt-caliber/first.pdb new file mode 100644 index 0000000000000000000000000000000000000000..71d596a288a82299863bab4ce30d4981a61cb850 --- /dev/null +++ b/regtest/isdb/rt-caliber/first.pdb @@ -0,0 +1,137 @@ +REMARK GENERATED BY TRJCONV +TITLE Macromolecule t= 0.00000 +MODEL 1 +ATOM 1 N GLY 1 8.280 2.110 -6.520 1.00 0.00 +ATOM 2 CA GLY 1 8.580 1.210 -7.670 1.00 0.00 +ATOM 3 C GLY 1 7.290 0.720 -8.320 1.00 0.00 +ATOM 4 O GLY 1 6.570 1.490 -8.930 1.00 0.00 +ATOM 5 N GLU 2 7.010 -0.550 -8.160 1.00 0.00 +ATOM 6 CA GLU 2 5.770 -1.110 -8.750 1.00 0.00 +ATOM 7 C GLU 2 4.690 -1.210 -7.670 1.00 0.00 +ATOM 8 O GLU 2 5.000 -1.440 -6.520 1.00 0.00 +ATOM 9 CB GLU 2 6.050 -2.490 -9.310 1.00 0.00 +ATOM 10 CG GLU 2 7.240 -2.410 -10.270 1.00 0.00 +ATOM 11 CD GLU 2 6.820 -2.930 -11.640 1.00 0.00 +ATOM 12 OE1 GLU 2 6.250 -2.140 -12.380 1.00 0.00 +ATOM 13 OE2 GLU 2 7.080 -4.100 -11.880 1.00 0.00 +ATOM 14 N TRP 3 3.460 -1.060 -8.070 1.00 0.00 +ATOM 15 CA TRP 3 2.340 -1.130 -7.080 1.00 0.00 +ATOM 16 C TRP 3 1.300 -2.160 -7.510 1.00 0.00 +ATOM 17 O TRP 3 0.940 -2.230 -8.660 1.00 0.00 +ATOM 18 CB TRP 3 1.700 0.240 -6.980 1.00 0.00 +ATOM 19 CG TRP 3 2.570 1.130 -6.090 1.00 0.00 +ATOM 20 CD1 TRP 3 3.740 1.620 -6.480 1.00 0.00 +ATOM 21 CD2 TRP 3 2.310 1.480 -4.830 1.00 0.00 +ATOM 22 CE2 TRP 3 3.340 2.260 -4.320 1.00 0.00 +ATOM 23 CE3 TRP 3 1.220 1.180 -4.010 1.00 0.00 +ATOM 24 NE1 TRP 3 4.180 2.300 -5.410 1.00 0.00 +ATOM 25 CZ2 TRP 3 3.280 2.720 -3.030 1.00 0.00 +ATOM 26 CZ3 TRP 3 1.180 1.650 -2.720 1.00 0.00 +ATOM 27 CH2 TRP 3 2.210 2.420 -2.220 1.00 0.00 +ATOM 28 N THR 4 0.850 -2.940 -6.560 1.00 0.00 +ATOM 29 CA THR 4 -0.180 -3.980 -6.870 1.00 0.00 +ATOM 30 C THR 4 -1.440 -3.780 -6.020 1.00 0.00 +ATOM 31 O THR 4 -1.430 -3.040 -5.060 1.00 0.00 +ATOM 32 CB THR 4 0.420 -5.360 -6.590 1.00 0.00 +ATOM 33 CG2 THR 4 1.440 -5.770 -7.650 1.00 0.00 +ATOM 34 OG1 THR 4 1.170 -5.200 -5.390 1.00 0.00 +ATOM 35 N TYR 5 -2.490 -4.450 -6.400 1.00 0.00 +ATOM 36 CA TYR 5 -3.770 -4.340 -5.650 1.00 0.00 +ATOM 37 C TYR 5 -4.430 -5.720 -5.570 1.00 0.00 +ATOM 38 O TYR 5 -4.610 -6.370 -6.580 1.00 0.00 +ATOM 39 CB TYR 5 -4.680 -3.350 -6.410 1.00 0.00 +ATOM 40 CG TYR 5 -5.830 -2.820 -5.510 1.00 0.00 +ATOM 41 CD1 TYR 5 -6.710 -3.670 -4.840 1.00 0.00 +ATOM 42 CD2 TYR 5 -6.010 -1.460 -5.380 1.00 0.00 +ATOM 43 CE1 TYR 5 -7.720 -3.150 -4.060 1.00 0.00 +ATOM 44 CE2 TYR 5 -7.020 -0.950 -4.590 1.00 0.00 +ATOM 45 CZ TYR 5 -7.880 -1.790 -3.930 1.00 0.00 +ATOM 46 OH TYR 5 -8.900 -1.280 -3.150 1.00 0.00 +ATOM 47 N ASP 6 -4.770 -6.130 -4.380 1.00 0.00 +ATOM 48 CA ASP 6 -5.430 -7.460 -4.210 1.00 0.00 +ATOM 49 C ASP 6 -6.920 -7.280 -3.900 1.00 0.00 +ATOM 50 O ASP 6 -7.350 -7.470 -2.780 1.00 0.00 +ATOM 51 CB ASP 6 -4.760 -8.210 -3.060 1.00 0.00 +ATOM 52 CG ASP 6 -5.080 -9.700 -3.170 1.00 0.00 +ATOM 53 OD1 ASP 6 -5.980 -10.010 -3.930 1.00 0.00 +ATOM 54 OD2 ASP 6 -4.400 -10.450 -2.480 1.00 0.00 +ATOM 55 N ASP 7 -7.680 -6.900 -4.890 1.00 0.00 +ATOM 56 CA ASP 7 -9.130 -6.690 -4.660 1.00 0.00 +ATOM 57 C ASP 7 -9.820 -8.020 -4.310 1.00 0.00 +ATOM 58 O ASP 7 -11.020 -8.090 -4.170 1.00 0.00 +ATOM 59 CB ASP 7 -9.770 -6.110 -5.920 1.00 0.00 +ATOM 60 CG ASP 7 -11.150 -5.540 -5.590 1.00 0.00 +ATOM 61 OD1 ASP 7 -11.200 -4.800 -4.620 1.00 0.00 +ATOM 62 OD2 ASP 7 -12.060 -5.890 -6.310 1.00 0.00 +ATOM 63 N ALA 8 -9.030 -9.050 -4.170 1.00 0.00 +ATOM 64 CA ALA 8 -9.620 -10.380 -3.830 1.00 0.00 +ATOM 65 C ALA 8 -9.730 -10.540 -2.310 1.00 0.00 +ATOM 66 O ALA 8 -10.510 -11.340 -1.830 1.00 0.00 +ATOM 67 CB ALA 8 -8.720 -11.490 -4.390 1.00 0.00 +ATOM 68 N THR 9 -8.950 -9.780 -1.590 1.00 0.00 +ATOM 69 CA THR 9 -9.010 -9.880 -0.100 1.00 0.00 +ATOM 70 C THR 9 -8.770 -8.500 0.530 1.00 0.00 +ATOM 71 O THR 9 -8.460 -8.400 1.700 1.00 0.00 +ATOM 72 CB THR 9 -7.930 -10.860 0.380 1.00 0.00 +ATOM 73 CG2 THR 9 -8.070 -12.230 -0.280 1.00 0.00 +ATOM 74 OG1 THR 9 -6.710 -10.320 -0.120 1.00 0.00 +ATOM 75 N LYS 10 -8.900 -7.480 -0.270 1.00 0.00 +ATOM 76 CA LYS 10 -8.680 -6.100 0.260 1.00 0.00 +ATOM 77 C LYS 10 -7.290 -5.990 0.900 1.00 0.00 +ATOM 78 O LYS 10 -7.160 -5.710 2.070 1.00 0.00 +ATOM 79 CB LYS 10 -9.750 -5.780 1.300 1.00 0.00 +ATOM 80 CG LYS 10 -11.140 -6.040 0.700 1.00 0.00 +ATOM 81 CD LYS 10 -12.200 -5.830 1.780 1.00 0.00 +ATOM 82 CE LYS 10 -12.720 -7.190 2.240 1.00 0.00 +ATOM 83 NZ LYS 10 -11.590 -8.060 2.650 1.00 0.00 +ATOM 84 N THR 11 -6.290 -6.240 0.100 1.00 0.00 +ATOM 85 CA THR 11 -4.890 -6.160 0.620 1.00 0.00 +ATOM 86 C THR 11 -3.990 -5.450 -0.410 1.00 0.00 +ATOM 87 O THR 11 -3.760 -5.970 -1.480 1.00 0.00 +ATOM 88 CB THR 11 -4.380 -7.600 0.830 1.00 0.00 +ATOM 89 CG2 THR 11 -3.310 -7.680 1.910 1.00 0.00 +ATOM 90 OG1 THR 11 -5.500 -8.310 1.350 1.00 0.00 +ATOM 91 N PHE 12 -3.510 -4.280 -0.060 1.00 0.00 +ATOM 92 CA PHE 12 -2.630 -3.560 -1.020 1.00 0.00 +ATOM 93 C PHE 12 -1.220 -4.150 -0.950 1.00 0.00 +ATOM 94 O PHE 12 -0.880 -4.820 0.000 1.00 0.00 +ATOM 95 CB PHE 12 -2.570 -2.090 -0.650 1.00 0.00 +ATOM 96 CG PHE 12 -3.240 -1.260 -1.730 1.00 0.00 +ATOM 97 CD1 PHE 12 -2.550 -0.890 -2.870 1.00 0.00 +ATOM 98 CD2 PHE 12 -4.550 -0.870 -1.590 1.00 0.00 +ATOM 99 CE1 PHE 12 -3.170 -0.130 -3.850 1.00 0.00 +ATOM 100 CE2 PHE 12 -5.160 -0.110 -2.550 1.00 0.00 +ATOM 101 CZ PHE 12 -4.480 0.270 -3.670 1.00 0.00 +ATOM 102 N THR 13 -0.430 -3.890 -1.950 1.00 0.00 +ATOM 103 CA THR 13 0.950 -4.430 -1.920 1.00 0.00 +ATOM 104 C THR 13 1.850 -3.650 -2.880 1.00 0.00 +ATOM 105 O THR 13 1.460 -3.320 -3.970 1.00 0.00 +ATOM 106 CB THR 13 0.920 -5.900 -2.350 1.00 0.00 +ATOM 107 CG2 THR 13 2.310 -6.520 -2.340 1.00 0.00 +ATOM 108 OG1 THR 13 0.200 -6.560 -1.320 1.00 0.00 +ATOM 109 N VAL 14 3.040 -3.360 -2.420 1.00 0.00 +ATOM 110 CA VAL 14 4.000 -2.610 -3.290 1.00 0.00 +ATOM 111 C VAL 14 5.370 -3.270 -3.190 1.00 0.00 +ATOM 112 O VAL 14 5.700 -3.860 -2.180 1.00 0.00 +ATOM 113 CB VAL 14 4.080 -1.150 -2.800 1.00 0.00 +ATOM 114 CG1 VAL 14 4.750 -1.120 -1.430 1.00 0.00 +ATOM 115 CG2 VAL 14 4.900 -0.320 -3.780 1.00 0.00 +ATOM 116 N THR 15 6.140 -3.170 -4.230 1.00 0.00 +ATOM 117 CA THR 15 7.490 -3.810 -4.220 1.00 0.00 +ATOM 118 C THR 15 8.500 -2.920 -4.940 1.00 0.00 +ATOM 119 O THR 15 8.410 -2.720 -6.130 1.00 0.00 +ATOM 120 CB THR 15 7.390 -5.160 -4.930 1.00 0.00 +ATOM 121 CG2 THR 15 8.580 -6.060 -4.620 1.00 0.00 +ATOM 122 OG1 THR 15 6.270 -5.810 -4.340 1.00 0.00 +ATOM 123 N GLU 16 9.450 -2.420 -4.200 1.00 0.00 +ATOM 124 CA GLU 16 10.480 -1.540 -4.840 1.00 0.00 +ATOM 125 C GLU 16 11.240 -2.320 -5.920 1.00 0.00 +ATOM 126 O GLU 16 12.090 -1.690 -6.540 1.00 0.00 +ATOM 127 CB GLU 16 11.470 -1.060 -3.780 1.00 0.00 +ATOM 128 CG GLU 16 11.710 0.440 -3.960 1.00 0.00 +ATOM 129 CD GLU 16 12.840 0.880 -3.020 1.00 0.00 +ATOM 130 OE1 GLU 16 12.820 0.420 -1.900 1.00 0.00 +ATOM 131 OE2 GLU 16 13.660 1.660 -3.490 1.00 0.00 +ATOM 132 OXT GLU 16 10.930 -3.490 -6.070 1.00 0.00 +TER +ENDMDL diff --git a/regtest/isdb/rt-caliber/plumed.dat b/regtest/isdb/rt-caliber/plumed.dat new file mode 100644 index 0000000000000000000000000000000000000000..4c8eb318d9240b912ab1600279cb9d58657283f3 --- /dev/null +++ b/regtest/isdb/rt-caliber/plumed.dat @@ -0,0 +1,93 @@ +MOLINFO STRUCTURE=first.pdb + +WHOLEMOLECULES ENTITY0=1-132 + +CONTACTMAP ... +ATOMS1=1,115 SWITCH1={Q REF=0.498366 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT1=0.013699 +ATOMS2=1,124 SWITCH2={Q REF=0.458093 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT2=0.013699 +ATOMS3=1,128 SWITCH3={Q REF=0.459428 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT3=0.013699 +ATOMS4=2,119 SWITCH4={Q REF=0.422438 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT4=0.013699 +ATOMS5=2,124 SWITCH5={Q REF=0.437966 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT5=0.013699 +ATOMS6=2,126 SWITCH6={Q REF=0.469116 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT6=0.013699 +ATOMS7=5,116 SWITCH7={Q REF=0.480273 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT7=0.013699 +ATOMS8=5,119 SWITCH8={Q REF=0.328478 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT8=0.013699 +ATOMS9=5,124 SWITCH9={Q REF=0.490341 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT9=0.013699 +ATOMS10=5,125 SWITCH10={Q REF=0.510327 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT10=0.013699 +ATOMS11=8,105 SWITCH11={Q REF=0.475063 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT11=0.013699 +ATOMS12=8,110 SWITCH12={Q REF=0.357796 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT12=0.013699 +ATOMS13=8,115 SWITCH13={Q REF=0.296176 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT13=0.013699 +ATOMS14=8,116 SWITCH14={Q REF=0.308814 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT14=0.013699 +ATOMS15=8,119 SWITCH15={Q REF=0.366314 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT15=0.013699 +ATOMS16=9,116 SWITCH16={Q REF=0.512610 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT16=0.013699 +ATOMS17=9,119 SWITCH17={Q REF=0.396672 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT17=0.013699 +ATOMS18=9,120 SWITCH18={Q REF=0.530178 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT18=0.013699 +ATOMS19=10,119 SWITCH19={Q REF=0.431330 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT19=0.013699 +ATOMS20=10,132 SWITCH20={Q REF=0.569408 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT20=0.013699 +ATOMS21=15,105 SWITCH21={Q REF=0.390418 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT21=0.013699 +ATOMS22=15,110 SWITCH22={Q REF=0.439433 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT22=0.013699 +ATOMS23=15,113 SWITCH23={Q REF=0.462022 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT23=0.013699 +ATOMS24=15,115 SWITCH24={Q REF=0.425437 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT24=0.013699 +ATOMS25=18,99 SWITCH25={Q REF=0.580092 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT25=0.013699 +ATOMS26=21,113 SWITCH26={Q REF=0.376440 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT26=0.013699 +ATOMS27=22,115 SWITCH27={Q REF=0.306294 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT27=0.013699 +ATOMS28=23,97 SWITCH28={Q REF=0.444943 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT28=0.013699 +ATOMS29=23,102 SWITCH29={Q REF=0.571586 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT29=0.013699 +ATOMS30=23,105 SWITCH30={Q REF=0.450657 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT30=0.013699 +ATOMS31=25,114 SWITCH31={Q REF=0.441209 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT31=0.013699 +ATOMS32=26,95 SWITCH32={Q REF=0.568639 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT32=0.013699 +ATOMS33=28,105 SWITCH33={Q REF=0.268786 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT33=0.013699 +ATOMS34=31,87 SWITCH34={Q REF=0.517979 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT34=0.013699 +ATOMS35=31,92 SWITCH35={Q REF=0.424641 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT35=0.013699 +ATOMS36=31,97 SWITCH36={Q REF=0.326696 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT36=0.013699 +ATOMS37=31,102 SWITCH37={Q REF=0.337559 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT37=0.013699 +ATOMS38=31,105 SWITCH38={Q REF=0.310139 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT38=0.013699 +ATOMS39=31,106 SWITCH39={Q REF=0.458761 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT39=0.013699 +ATOMS40=33,122 SWITCH40={Q REF=0.585548 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT40=0.013699 +ATOMS41=34,103 SWITCH41={Q REF=0.356121 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT41=0.013699 +ATOMS42=34,105 SWITCH42={Q REF=0.237379 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT42=0.013699 +ATOMS43=34,106 SWITCH43={Q REF=0.312955 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT43=0.013699 +ATOMS44=34,107 SWITCH44={Q REF=0.351347 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT44=0.013699 +ATOMS45=34,122 SWITCH45={Q REF=0.524258 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT45=0.013699 +ATOMS46=36,87 SWITCH46={Q REF=0.447726 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT46=0.013699 +ATOMS47=36,92 SWITCH47={Q REF=0.483166 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT47=0.013699 +ATOMS48=36,97 SWITCH48={Q REF=0.459557 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT48=0.013699 +ATOMS49=36,102 SWITCH49={Q REF=0.500481 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT49=0.013699 +ATOMS50=39,99 SWITCH50={Q REF=0.438202 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT50=0.013699 +ATOMS51=40,92 SWITCH51={Q REF=0.556307 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT51=0.013699 +ATOMS52=41,87 SWITCH52={Q REF=0.502813 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT52=0.013699 +ATOMS53=42,101 SWITCH53={Q REF=0.287366 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT53=0.013699 +ATOMS54=43,76 SWITCH54={Q REF=0.531851 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT54=0.013699 +ATOMS55=43,84 SWITCH55={Q REF=0.537574 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT55=0.013699 +ATOMS56=44,100 SWITCH56={Q REF=0.288562 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT56=0.013699 +ATOMS57=46,76 SWITCH57={Q REF=0.590838 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT57=0.013699 +ATOMS58=47,84 SWITCH58={Q REF=0.473211 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT58=0.013699 +ATOMS59=47,87 SWITCH59={Q REF=0.307501 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT59=0.013699 +ATOMS60=47,92 SWITCH60={Q REF=0.474069 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT60=0.013699 +ATOMS61=47,96 SWITCH61={Q REF=0.575155 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT61=0.013699 +ATOMS62=47,108 SWITCH62={Q REF=0.585230 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT62=0.013699 +ATOMS63=50,75 SWITCH63={Q REF=0.295003 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT63=0.013699 +ATOMS64=50,84 SWITCH64={Q REF=0.330619 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT64=0.013699 +ATOMS65=50,87 SWITCH65={Q REF=0.410221 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT65=0.013699 +ATOMS66=50,90 SWITCH66={Q REF=0.460272 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT66=0.013699 +ATOMS67=51,84 SWITCH67={Q REF=0.402584 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT67=0.013699 +ATOMS68=51,85 SWITCH68={Q REF=0.421448 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT68=0.013699 +ATOMS69=51,87 SWITCH69={Q REF=0.291788 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT69=0.013699 +ATOMS70=51,88 SWITCH70={Q REF=0.395583 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT70=0.013699 +ATOMS71=51,108 SWITCH71={Q REF=0.550924 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT71=0.013699 +ATOMS72=52,75 SWITCH72={Q REF=0.528496 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT72=0.013699 +ATOMS73=54,88 SWITCH73={Q REF=0.436795 BETA=50.0 LAMBDA=1.8 R_0=0.01} WEIGHT73=0.013699 +SUM +LABEL=pairs +... CONTACTMAP + +CALIBER ... + ARG=pairs + KAPPA=1600000 + FILE=Qaverages.dat + LABEL=cal +... CALIBER + + +PRINT ARG=pairs,cal.x0_0,cal.bias FILE=COLVAR STRIDE=1 + +ENDPLUMED diff --git a/regtest/isdb/rt-caliber/trajout.0.xtc b/regtest/isdb/rt-caliber/trajout.0.xtc new file mode 100644 index 0000000000000000000000000000000000000000..be6c5a92a10cfecea8403fd2047d3851806afa26 Binary files /dev/null and b/regtest/isdb/rt-caliber/trajout.0.xtc differ diff --git a/regtest/isdb/rt-caliber/trajout.1.xtc b/regtest/isdb/rt-caliber/trajout.1.xtc new file mode 100644 index 0000000000000000000000000000000000000000..a2380fe8aaefcd7db8894dd0f5a04251df5ab0c1 Binary files /dev/null and b/regtest/isdb/rt-caliber/trajout.1.xtc differ diff --git a/regtest/isdb/rt-caliber/trajout.2.xtc b/regtest/isdb/rt-caliber/trajout.2.xtc new file mode 100644 index 0000000000000000000000000000000000000000..0b2c4667c25361bf9bdbd19208109aab9be21a54 Binary files /dev/null and b/regtest/isdb/rt-caliber/trajout.2.xtc differ diff --git a/regtest/isdb/rt-caliber/trajout.3.xtc b/regtest/isdb/rt-caliber/trajout.3.xtc new file mode 100644 index 0000000000000000000000000000000000000000..fadfb72b534d98f5f45318616b5fe12206812d4c Binary files /dev/null and b/regtest/isdb/rt-caliber/trajout.3.xtc differ diff --git a/src/isdb/Caliber.cpp b/src/isdb/Caliber.cpp new file mode 100644 index 0000000000000000000000000000000000000000..5d90860d24073bcc6d62c352af34986b4297f5e4 --- /dev/null +++ b/src/isdb/Caliber.cpp @@ -0,0 +1,344 @@ +/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + Copyright (c) 2017,2018 The plumed team + (see the PEOPLE file at the root of the distribution for a list of names) + + See http://www.plumed.org for more information. + + This file is part of plumed, version 2. + + plumed is free software: you can redistribute it and/or modify + it under the terms of the GNU Lesser General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + plumed is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public License + along with plumed. If not, see <http://www.gnu.org/licenses/>. ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */ +#include "bias/Bias.h" +#include "bias/ActionRegister.h" +#include "core/Atoms.h" +#include "core/PlumedMain.h" +#include <fstream> + +using namespace std; + +namespace PLMD { +namespace isdb { + +//+PLUMEDOC BIAS CALIBER +/* +Add a time-dependent, harmonic restraint on one or more variables. +This allows implementing a maximum caliber restraint on one or more experimental time serie by replica-averaged restrained simulations. +See \cite{}. + +The time resolved experiments are read from a text file and intermediate values are obtained by splines. + +\par Examples + +In the following example a restraint is applied on the time evolution of a saxs spectrum + +MOLINFO STRUCTURE=first.pdb + +# Define saxs variable +SAXS ... +LABEL=saxs +ATOMISTIC +ATOMS=1-436 +QVALUE1=0.02 # Q-value at which calculate the scattering +QVALUE2=0.0808 +QVALUE3=0.1264 +QVALUE4=0.1568 +QVALUE5=0.172 +QVALUE6=0.1872 +QVALUE7=0.2176 +QVALUE8=0.2328 +QVALUE9=0.248 +QVALUE10=0.2632 +QVALUE11=0.2936 +QVALUE12=0.3088 +QVALUE13=0.324 +QVALUE14=0.3544 +QVALUE15=0.4 +... SAXS + + +#define the caliber restraint +CALIBER ... + ARG=(saxs\.q_.*) + FILE=expsaxs.dat + KAPPA=10 + LABEL=cal0 + STRIDE=10 + REGRES_ZERO=200 + AVERAGING=200 +... CALIBER + +In particular the file expsaxs.dat contains the time traces for the 15 intensities at the selected scattering lengths, organised as time, q_1, etc. +The strenght of the bias is automatically evaluated from the standard error of the mean over AVERAGING steps and multiplied by KAPPA. This is usefull when working with multiple experimental data +Because SAXS is usually defined irrespectively of a scaling factor the scaling is evaluated from a linear fit every REGRES_ZERO step. Alternatively it can be given as a fixed constant as SCALE. +The bias is here applied every 10th steps. + +*/ +//+ENDPLUMEDOC + + +class Caliber : public bias::Bias { +public: + explicit Caliber(const ActionOptions&); + void calculate(); + static void registerKeywords( Keywords& keys ); +private: + vector<double> time; + vector< vector<double> > var; + vector< vector<double> > dvar; + double mult; + double scale_; + bool master; + unsigned replica_; + unsigned nrep_; + // scale and offset regression + bool doregres_zero_; + int nregres_zero_; + // force constant + unsigned optsigmamean_stride_; + vector<double> sigma_mean2_; + vector< vector<double> > sigma_mean2_last_; + vector<Value*> x0comp; + vector<Value*> kcomp; + vector<Value*> mcomp; + Value* valueScale; + + void get_sigma_mean(const double fact, const vector<double> &mean); + void replica_averaging(const double fact, vector<double> &mean); + double getSpline(const unsigned iarg); + void do_regression_zero(const vector<double> &mean); +}; + +PLUMED_REGISTER_ACTION(Caliber,"CALIBER") + +void Caliber::registerKeywords( Keywords& keys ) { + Bias::registerKeywords(keys); + keys.use("ARG"); + keys.addFlag("NOENSEMBLE",false,"don't perform any replica-averaging"); + keys.add("compulsory","FILE","the name of the file containing the time-resolved values"); + keys.add("compulsory","KAPPA","a force constant, this can be use to scale a constant estimanted on-the-fly using AVERAGING"); + keys.add("optional","AVERAGING", "Stride for calculation of the optimum kappa, if 0 only KAPPA is used."); + keys.add("compulsory","TSCALE","1.0","Apply a time scaling on the experimental time scale"); + keys.add("compulsory","SCALE","1.0","Apply a constant scaling on the data provided as arguments"); + keys.add("optional","REGRES_ZERO","stride for regression with zero offset"); + keys.addOutputComponent("x0","default","the instantaneous value of the center of the potential"); + keys.addOutputComponent("mean","default","the current average value of the calculated observable"); + keys.addOutputComponent("kappa","default","the current force constant"); + keys.addOutputComponent("scale","REGRES_ZERO","the current scaling constant"); +} + +Caliber::Caliber(const ActionOptions&ao): + PLUMED_BIAS_INIT(ao), + mult(0), + scale_(1), + doregres_zero_(false), + nregres_zero_(0), + optsigmamean_stride_(0) +{ + parse("KAPPA",mult); + string filename; + parse("FILE",filename); + if( filename.length()==0 ) error("No external variable file was specified"); + unsigned averaging=0; + parse("AVERAGING", averaging); + if(averaging>0) optsigmamean_stride_ = averaging; + double tscale=1.0; + parse("TSCALE", tscale); + if(tscale<=0.) error("The time scale factor must be greater than 0."); + parse("SCALE", scale_); + if(scale_==0.) error("The time scale factor cannot be 0."); + // regression with zero intercept + parse("REGRES_ZERO", nregres_zero_); + if(nregres_zero_>0) { + // set flag + doregres_zero_=true; + log.printf(" doing regression with zero intercept with stride: %d\n", nregres_zero_); + } + + + bool noensemble = false; + parseFlag("NOENSEMBLE", noensemble); + + checkRead(); + + // set up replica stuff + master = (comm.Get_rank()==0); + if(master) { + nrep_ = multi_sim_comm.Get_size(); + replica_ = multi_sim_comm.Get_rank(); + if(noensemble) nrep_ = 1; + } else { + nrep_ = 0; + replica_ = 0; + } + comm.Sum(&nrep_,1); + comm.Sum(&replica_,1); + + const unsigned narg = getNumberOfArguments(); + sigma_mean2_.resize(narg,1); + sigma_mean2_last_.resize(narg); + for(unsigned j=0; j<narg; j++) sigma_mean2_last_[j].push_back(0.000001); + + log.printf(" Time resolved data from file %s\n",filename.c_str()); + std::ifstream varfile(filename.c_str()); + if(varfile.fail()) error("Cannot open "+filename); + var.resize(narg); + dvar.resize(narg); + while (!varfile.eof()) { + double tempT, tempVar; + varfile >> tempT; + time.push_back(tempT/tscale); + for(unsigned i=0; i<narg; i++) { + varfile >> tempVar; + var[i].push_back(tempVar); + } + } + varfile.close(); + + const double deltat = time[1] - time[0]; + for(unsigned i=0; i<narg; i++) { + for(unsigned j=0; j<var[i].size(); j++) { + if(j==0) dvar[i].push_back((var[i][j+1] - var[i][j])/(deltat)); + else if(j==var[i].size()-1) dvar[i].push_back((var[i][j] - var[i][j-1])/(deltat)); + else dvar[i].push_back((var[i][j+1] - var[i][j-1])/(2.*deltat)); + } + } + + for(unsigned i=0; i<narg; i++) { + std::string num; Tools::convert(i,num); + addComponent("x0_"+num); componentIsNotPeriodic("x0_"+num); x0comp.push_back(getPntrToComponent("x0_"+num)); + addComponent("kappa_"+num); componentIsNotPeriodic("kappa_"+num); kcomp.push_back(getPntrToComponent("kappa_"+num)); + addComponent("mean_"+num); componentIsNotPeriodic("mean_"+num); mcomp.push_back(getPntrToComponent("mean_"+num)); + } + + if(doregres_zero_) { + addComponent("scale"); + componentIsNotPeriodic("scale"); + valueScale=getPntrToComponent("scale"); + } +} + +void Caliber::get_sigma_mean(const double fact, const vector<double> &mean) +{ + const unsigned narg = getNumberOfArguments(); + const double dnrep = static_cast<double>(nrep_); + + if(sigma_mean2_last_[0].size()==optsigmamean_stride_) for(unsigned i=0; i<narg; ++i) sigma_mean2_last_[i].erase(sigma_mean2_last_[i].begin()); + vector<double> sigma_mean2_now(narg,0); + if(master) { + for(unsigned i=0; i<narg; ++i) { + double tmp = getArgument(i)-mean[i]; + sigma_mean2_now[i] = fact*tmp*tmp; + } + if(nrep_>1) multi_sim_comm.Sum(&sigma_mean2_now[0], narg); + } + comm.Sum(&sigma_mean2_now[0], narg); + + for(unsigned i=0; i<narg; ++i) { + sigma_mean2_last_[i].push_back(sigma_mean2_now[i]/dnrep); + sigma_mean2_[i] = *max_element(sigma_mean2_last_[i].begin(), sigma_mean2_last_[i].end()); + } +} + +void Caliber::replica_averaging(const double fact, vector<double> &mean) +{ + const unsigned narg = getNumberOfArguments(); + if(master) { + for(unsigned i=0; i<narg; ++i) mean[i] = fact*getArgument(i); + if(nrep_>1) multi_sim_comm.Sum(&mean[0], narg); + } + comm.Sum(&mean[0], narg); +} + +double Caliber::getSpline(const unsigned iarg) +{ + const double deltat = time[1] - time[0]; + const int tindex = static_cast<int>(getTime()/deltat); + + unsigned start, end; + start=tindex; + if(tindex+1<var[iarg].size()) end=tindex+2; + else end=var[iarg].size(); + + double value=0; + for(unsigned ipoint=start; ipoint<end; ++ipoint) { + double grid=var[iarg][ipoint]; + double dder=dvar[iarg][ipoint]; + double yy=0.; + if(fabs(grid)>0.0000001) yy=-dder/grid; + + int x0=1; + if(ipoint==tindex) x0=0; + + double X=fabs((getTime()-time[tindex])/deltat-(double)x0); + double X2=X*X; + double X3=X2*X; + double C=(1.0-3.0*X2+2.0*X3) - (x0?-1.0:1.0)*yy*(X-2.0*X2+X3)*deltat; + + value+=grid*C; + } + return value; +} + +void Caliber::do_regression_zero(const vector<double> &mean) +{ +// parameters[i] = scale_ * mean[i]: find scale_ with linear regression + double num = 0.0; + double den = 0.0; + for(unsigned i=0; i<getNumberOfArguments(); ++i) { + num += mean[i] * getSpline(i); + den += mean[i] * mean[i]; + } + if(den>0) { + scale_ = num / den; + } else { + scale_ = 1.0; + } +} + +void Caliber::calculate() +{ + const unsigned narg = getNumberOfArguments(); + const double dnrep = static_cast<double>(nrep_); + const double fact = 1.0/dnrep; + + vector<double> mean(narg,0); + vector<double> dmean_x(narg,fact); + replica_averaging(fact, mean); + if(optsigmamean_stride_>0) get_sigma_mean(fact, mean); + + // in case of regression with zero intercept, calculate scale + if(doregres_zero_ && getStep()%nregres_zero_==0) do_regression_zero(mean); + + double ene=0; + for(unsigned i=0; i<narg; ++i) { + const double x0 = getSpline(i); + const double kappa = mult*dnrep/sigma_mean2_[i]; + const double cv=difference(i,x0,scale_*mean[i]); + const double f=-kappa*cv*dmean_x[i]/scale_; + setOutputForce(i,f); + ene+=0.5*kappa*cv*cv; + x0comp[i]->set(x0); + kcomp[i]->set(kappa); + mcomp[i]->set(mean[i]); + } + + if(doregres_zero_) valueScale->set(scale_); + + setBias(ene); +} + +} +} + + diff --git a/src/isdb/Metainference.cpp b/src/isdb/Metainference.cpp index b46589ca2ef40a59030cfe5968b4bc7b76359e2e..ad9c8c55570f72e6d8d7b3de8d954bc61b901466 100644 --- a/src/isdb/Metainference.cpp +++ b/src/isdb/Metainference.cpp @@ -173,6 +173,9 @@ class Metainference : public bias::Bias double offset_min_; double offset_max_; double Doffset_; + // scale and offset regression + bool doregres_zero_; + int nregres_zero_; // sigma is data uncertainty vector<double> sigma_; vector<double> sigma_min_; @@ -251,7 +254,8 @@ class Metainference : public bias::Bias void get_weights(const unsigned iselect, double &fact, double &var_fact); void replica_averaging(const double fact, std::vector<double> &mean, std::vector<double> &dmean_b); void get_sigma_mean(const unsigned iselect, const double fact, const double var_fact, const vector<double> &mean); - void writeStatus(); + void writeStatus(); + void do_regression_zero(const vector<double> &mean); public: explicit Metainference(const ActionOptions&); @@ -287,6 +291,7 @@ void Metainference::registerKeywords(Keywords& keys) { keys.add("optional","OFFSET_MIN","minimum value of the offset"); keys.add("optional","OFFSET_MAX","maximum value of the offset"); keys.add("optional","DOFFSET","maximum MC move of the offset"); + keys.add("optional","REGRES_ZERO","stride for regression with zero offset"); keys.add("compulsory","SIGMA0","1.0","initial value of the uncertainty parameter"); keys.add("compulsory","SIGMA_MIN","0.0","minimum value of the uncertainty parameter"); keys.add("compulsory","SIGMA_MAX","10.","maximum value of the uncertainty parameter"); @@ -309,8 +314,8 @@ void Metainference::registerKeywords(Keywords& keys) { keys.addOutputComponent("acceptScale", "SCALEDATA", "MC acceptance"); keys.addOutputComponent("weight", "REWEIGHT", "weights of the weighted average"); keys.addOutputComponent("biasDer", "REWEIGHT", "derivatives wrt the bias"); - keys.addOutputComponent("scale", "SCALEDATA", "scale parameter"); - keys.addOutputComponent("offset", "ADDOFFSET", "offset parameter"); + keys.addOutputComponent("scale", "default", "scale parameter"); + keys.addOutputComponent("offset", "default", "offset parameter"); keys.addOutputComponent("ftilde", "GENERIC", "ensemble average estimator"); } @@ -328,6 +333,8 @@ Metainference::Metainference(const ActionOptions&ao): offset_min_(1), offset_max_(-1), Doffset_(-1), + doregres_zero_(false), + nregres_zero_(0), Dftilde_(0.1), random(3), MCsteps_(1), @@ -506,6 +513,16 @@ Metainference::Metainference(const ActionOptions&ao): } } + // regression with zero intercept + parse("REGRES_ZERO", nregres_zero_); + if(nregres_zero_>0) { + // set flag + doregres_zero_=true; + // check if already sampling scale and offset + if(doscale_) error("REGRES_ZERO and SCALEDATA are mutually exclusive"); + if(dooffset_) error("REGRES_ZERO and ADDOFFSET are mutually exclusive"); + } + vector<double> readsigma; parseVector("SIGMA0",readsigma); if((noise_type_!=MGAUSS&&noise_type_!=MOUTLIERS&&noise_type_!=GENERIC)&&readsigma.size()>1) @@ -658,6 +675,9 @@ Metainference::Metainference(const ActionOptions&ao): } } + if(doregres_zero_) + log.printf(" doing regression with zero intercept with stride: %d\n", nregres_zero_); + log.printf(" number of experimental data points %u\n",narg); log.printf(" number of replicas %u\n",nrep_); log.printf(" initial data uncertainties"); @@ -680,7 +700,7 @@ Metainference::Metainference(const ActionOptions&ao): componentIsNotPeriodic("weight"); } - if(doscale_) { + if(doscale_ || doregres_zero_) { addComponent("scale"); componentIsNotPeriodic("scale"); valueScale=getPntrToComponent("scale"); @@ -788,7 +808,7 @@ double Metainference::getEnergySP(const vector<double> &mean, const vector<doubl } // add one single Jeffrey's prior and one normalisation per data point ene += 0.5*std::log(sss) + static_cast<double>(narg)*0.5*std::log(0.5*M_PI*M_PI/ss2); - if(doscale_) ene += 0.5*std::log(sss); + if(doscale_ || doregres_zero_) ene += 0.5*std::log(sss); if(dooffset_) ene += 0.5*std::log(sss); return kbt_ * ene; } @@ -808,7 +828,7 @@ double Metainference::getEnergySPE(const vector<double> &mean, const vector<doub const double dev = scale*mean[i]-parameters[i]+offset; const double a2 = 0.5*dev*dev + ss2; ene += 0.5*std::log(sss) + 0.5*std::log(0.5*M_PI*M_PI/ss2) + std::log(2.0*a2/(1.0-exp(-a2/sm2))); - if(doscale_) ene += 0.5*std::log(sss); + if(doscale_ || doregres_zero_) ene += 0.5*std::log(sss); if(dooffset_) ene += 0.5*std::log(sss); } } @@ -836,7 +856,7 @@ double Metainference::getEnergyMIGEN(const vector<double> &mean, const vector<do const double normm = -0.5*std::log(0.5/M_PI*inv_sm2); const double jeffreys = -0.5*std::log(2.*inv_sb2); ene += 0.5*devb*devb*inv_sb2 + 0.5*devm*devm*inv_sm2 + normb + normm + jeffreys; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; } } @@ -863,7 +883,7 @@ double Metainference::getEnergyGJ(const vector<double> &mean, const vector<doubl const double jeffreys = -0.5*std::log(2.*inv_sss); // add Jeffrey's prior in case one sigma for all data points + one normalisation per datapoint ene += jeffreys + static_cast<double>(narg)*normalisation; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; return kbt_ * ene; @@ -886,7 +906,7 @@ double Metainference::getEnergyGJE(const vector<double> &mean, const vector<doub const double normalisation = -0.5*std::log(0.5/M_PI*inv_s2); const double jeffreys = -0.5*std::log(2.*inv_sss); ene += 0.5*dev*dev*inv_s2 + normalisation + jeffreys; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; } } @@ -1133,7 +1153,7 @@ void Metainference::doMonteCarlo(const vector<double> &mean_) /* save the result of the sampling */ double accept = static_cast<double>(MCaccept_) / static_cast<double>(MCtrial_); valueAccept->set(accept); - if(doscale_) valueScale->set(scale_); + if(doscale_ || doregres_zero_) valueScale->set(scale_); if(dooffset_) valueOffset->set(offset_); if(doscale_||dooffset_) { accept = static_cast<double>(MCacceptScale_) / static_cast<double>(MCtrial_); @@ -1499,17 +1519,35 @@ void Metainference::replica_averaging(const double fact, vector<double> &mean, v if(firstTime) {ftilde_ = mean; firstTime = false;} } +void Metainference::do_regression_zero(const vector<double> &mean) +{ +// parameters[i] = scale_ * mean[i]: find scale_ with linear regression + double num = 0.0; + double den = 0.0; + for(unsigned i=0; i<parameters.size(); ++i) { + num += mean[i] * parameters[i]; + den += mean[i] * mean[i]; + } + if(den>0) { + scale_ = num / den; + } else { + scale_ = 1.0; + } +} + void Metainference::calculate() { + // get step + const long int step = getStep(); + unsigned iselect = 0; // set the value of selector for REM-like stuff if(selector_.length()>0) iselect = static_cast<unsigned>(plumed.passMap[selector_]); double fact = 0.0; double var_fact = 0.0; - + // get weights for ensemble average get_weights(iselect, fact, var_fact); - // calculate the mean vector<double> mean(narg,0); // this is the derivative of the mean with respect to the argument @@ -1518,12 +1556,12 @@ void Metainference::calculate() vector<double> dmean_b(narg,0); // calculate it replica_averaging(fact, mean, dmean_b); - + // calculate sigma mean get_sigma_mean(iselect, fact, var_fact, mean); - + // in case of regression with zero intercept, calculate scale + if(doregres_zero_ && step%nregres_zero_==0) do_regression_zero(mean); /* MONTE CARLO */ - const long int step = getStep(); if(step%MCstride_==0&&!getExchangeStep()) doMonteCarlo(mean); // calculate bias and forces diff --git a/src/isdb/MetainferenceBase.cpp b/src/isdb/MetainferenceBase.cpp index c94c7b017939794f1e0e38b8e64debfa92a0245a..737bd74fc6cee82fad30110f52e84cfa8ef70a9a 100644 --- a/src/isdb/MetainferenceBase.cpp +++ b/src/isdb/MetainferenceBase.cpp @@ -60,6 +60,7 @@ void MetainferenceBase::registerKeywords( Keywords& keys ) { keys.add("optional","OFFSET_MIN","minimum value of the offset"); keys.add("optional","OFFSET_MAX","maximum value of the offset"); keys.add("optional","DOFFSET","maximum MC move of the offset"); + keys.add("optional","REGRES_ZERO","stride for regression with zero offset"); keys.add("compulsory","SIGMA0","1.0","initial value of the uncertainty parameter"); keys.add("compulsory","SIGMA_MIN","0.0","minimum value of the uncertainty parameter"); keys.add("compulsory","SIGMA_MAX","10.","maximum value of the uncertainty parameter"); @@ -107,6 +108,8 @@ MetainferenceBase::MetainferenceBase(const ActionOptions&ao): offset_min_(1), offset_max_(-1), Doffset_(-1), + doregres_zero_(false), + nregres_zero_(0), Dftilde_(0.1), random(3), MCsteps_(1), @@ -255,6 +258,16 @@ MetainferenceBase::MetainferenceBase(const ActionOptions&ao): } } + // regression with zero intercept + parse("REGRES_ZERO", nregres_zero_); + if(nregres_zero_>0) { + // set flag + doregres_zero_=true; + // check if already sampling scale and offset + if(doscale_) error("REGRES_ZERO and SCALEDATA are mutually exclusive"); + if(dooffset_) error("REGRES_ZERO and ADDOFFSET are mutually exclusive"); + } + vector<double> readsigma; parseVector("SIGMA0",readsigma); if((noise_type_!=MGAUSS&&noise_type_!=MOUTLIERS&&noise_type_!=GENERIC)&&readsigma.size()>1) @@ -471,7 +484,7 @@ void MetainferenceBase::Initialise(const unsigned input) componentIsNotPeriodic("weight"); } - if(doscale_) { + if(doscale_ || doregres_zero_) { addComponent("scale"); componentIsNotPeriodic("scale"); valueScale=getPntrToComponent("scale"); @@ -626,7 +639,7 @@ double MetainferenceBase::getEnergySP(const vector<double> &mean, const vector<d } // add one single Jeffrey's prior and one normalisation per data point ene += 0.5*std::log(sss) + static_cast<double>(narg)*0.5*std::log(0.5*M_PI*M_PI/ss2); - if(doscale_) ene += 0.5*std::log(sss); + if(doscale_ || doregres_zero_) ene += 0.5*std::log(sss); if(dooffset_) ene += 0.5*std::log(sss); return kbt_ * ene; } @@ -646,7 +659,7 @@ double MetainferenceBase::getEnergySPE(const vector<double> &mean, const vector< const double dev = scale*mean[i]-parameters[i]+offset; const double a2 = 0.5*dev*dev + ss2; ene += 0.5*std::log(sss) + 0.5*std::log(0.5*M_PI*M_PI/ss2) + std::log(2.0*a2/(1.0-exp(-a2/sm2))); - if(doscale_) ene += 0.5*std::log(sss); + if(doscale_ || doregres_zero_) ene += 0.5*std::log(sss); if(dooffset_) ene += 0.5*std::log(sss); } } @@ -674,7 +687,7 @@ double MetainferenceBase::getEnergyMIGEN(const vector<double> &mean, const vecto const double normm = -0.5*std::log(0.5/M_PI*inv_sm2); const double jeffreys = -0.5*std::log(2.*inv_sb2); ene += 0.5*devb*devb*inv_sb2 + 0.5*devm*devm*inv_sm2 + normb + normm + jeffreys; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; } } @@ -701,7 +714,7 @@ double MetainferenceBase::getEnergyGJ(const vector<double> &mean, const vector<d const double jeffreys = -0.5*std::log(2.*inv_sss); // add Jeffrey's prior in case one sigma for all data points + one normalisation per datapoint ene += jeffreys + static_cast<double>(narg)*normalisation; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; return kbt_ * ene; @@ -724,7 +737,7 @@ double MetainferenceBase::getEnergyGJE(const vector<double> &mean, const vector< const double normalisation = -0.5*std::log(0.5/M_PI*inv_s2); const double jeffreys = -0.5*std::log(2.*inv_sss); ene += 0.5*dev*dev*inv_s2 + normalisation + jeffreys; - if(doscale_) ene += jeffreys; + if(doscale_ || doregres_zero_) ene += jeffreys; if(dooffset_) ene += jeffreys; } } @@ -973,7 +986,7 @@ void MetainferenceBase::doMonteCarlo(const vector<double> &mean_) /* save the result of the sampling */ double accept = static_cast<double>(MCaccept_) / static_cast<double>(MCtrial_); valueAccept->set(accept); - if(doscale_) valueScale->set(scale_); + if(doscale_ || doregres_zero_) valueScale->set(scale_); if(dooffset_) valueOffset->set(offset_); if(doscale_||dooffset_) { accept = static_cast<double>(MCacceptScale_) / static_cast<double>(MCtrial_); @@ -1338,6 +1351,22 @@ void MetainferenceBase::replica_averaging(const double fact, vector<double> &mea if(firstTime) {ftilde_ = mean; firstTime = false;} } +void MetainferenceBase::do_regression_zero(const vector<double> &mean) +{ +// parameters[i] = scale_ * mean[i]: find scale_ with linear regression + double num = 0.0; + double den = 0.0; + for(unsigned i=0; i<parameters.size(); ++i) { + num += mean[i] * parameters[i]; + den += mean[i] * mean[i]; + } + if(den>0) { + scale_ = num / den; + } else { + scale_ = 1.0; + } +} + double MetainferenceBase::getScore() { /* Metainference */ @@ -1358,6 +1387,9 @@ double MetainferenceBase::getScore() /* 3) calculates parameters */ get_sigma_mean(fact, var_fact, mean); + // in case of regression with zero intercept, calculate scale + if(doregres_zero_ && getStep()%nregres_zero_==0) do_regression_zero(mean); + /* 4) run monte carlo */ doMonteCarlo(mean); diff --git a/src/isdb/MetainferenceBase.h b/src/isdb/MetainferenceBase.h index b5b651d292186fcc5f657f64adcb98b7855575f1..aa3663cb1f7e1346916808886996016e8fd1bb97 100644 --- a/src/isdb/MetainferenceBase.h +++ b/src/isdb/MetainferenceBase.h @@ -83,6 +83,9 @@ private: double offset_min_; double offset_max_; double Doffset_; + // scale and offset regression + bool doregres_zero_; + int nregres_zero_; // sigma is data uncertainty std::vector<double> sigma_; std::vector<double> sigma_min_; @@ -164,6 +167,7 @@ private: void get_weights(double &fact, double &var_fact); void replica_averaging(const double fact, std::vector<double> &mean, std::vector<double> &dmean_b); void get_sigma_mean(const double fact, const double var_fact, const std::vector<double> &mean); + void do_regression_zero(const std::vector<double> &mean); void doMonteCarlo(const std::vector<double> &mean);