From bd67911a8d5a1d94099bbbb49485076154fd7b0c Mon Sep 17 00:00:00 2001
From: Gareth Tribello <gareth.tribello@gmail.com>
Date: Mon, 25 Mar 2013 20:28:30 +0000
Subject: [PATCH] Changes to make multicolvar faster when dealing with
coordination number like variables
The code now only adds derivatives for the quantities being accumulated
when the values are greater than TOL. In addition, all derivative
accumulation can now be done in vesselbase. (MergeDerivatives is still in
MultiColvar but only because this makes it faster.) This makes it easier to add new
object that inherit from ActionWithVessels. Lastly I split the two
parts of apply in MultiColvar to remove code duplication.
---
regtest/rt22/derivatives2.reference | 392 ++++++++++----------
regtest/rt22/plumed.dat | 4 +-
src/core/ActionAtomistic.cpp | 18 +
src/core/ActionAtomistic.h | 2 +
src/multicolvar/ActionVolume.cpp | 17 +-
src/multicolvar/ActionVolume.h | 7 -
src/multicolvar/MultiColvar.cpp | 139 ++++---
src/multicolvar/MultiColvar.h | 127 +++----
src/multicolvar/MultiColvarFunction.cpp | 17 +-
src/multicolvar/MultiColvarFunction.h | 11 +-
src/multicolvar/StoreCentralAtomsVessel.cpp | 4 +-
src/vesselbase/ActionWithVessel.cpp | 82 +++-
src/vesselbase/ActionWithVessel.h | 44 ++-
src/vesselbase/BridgeVessel.cpp | 7 +-
src/vesselbase/Moments.cpp | 4 +-
src/vesselbase/StoreValuesVessel.cpp | 35 +-
src/vesselbase/StoreValuesVessel.h | 10 +-
src/vesselbase/Vessel.cpp | 2 +-
18 files changed, 505 insertions(+), 417 deletions(-)
diff --git a/regtest/rt22/derivatives2.reference b/regtest/rt22/derivatives2.reference
index a30da2409..33e68b8bc 100644
--- a/regtest/rt22/derivatives2.reference
+++ b/regtest/rt22/derivatives2.reference
@@ -1,196 +1,196 @@
-#! FIELDS time parameter s2.mean s2n.mean
- 0.000000 0 0.3185 0.3185
- 0.000000 1 0.2157 0.2157
- 0.000000 2 0.1024 0.1024
- 0.000000 3 -0.3069 -0.3069
- 0.000000 4 0.4357 0.4357
- 0.000000 5 0.0314 0.0314
- 0.000000 6 -0.3808 -0.3808
- 0.000000 7 -0.3976 -0.3976
- 0.000000 8 0.2787 0.2787
- 0.000000 9 0.5133 0.5133
- 0.000000 10 -0.4440 -0.4440
- 0.000000 11 0.0141 0.0141
- 0.000000 12 0.3969 0.3969
- 0.000000 13 0.5493 0.5493
- 0.000000 14 0.0211 0.0211
- 0.000000 15 -0.3772 -0.3772
- 0.000000 16 0.2829 0.2829
- 0.000000 17 -0.2163 -0.2163
- 0.000000 18 -0.5101 -0.5101
- 0.000000 19 -0.4685 -0.4685
- 0.000000 20 -0.0317 -0.0317
- 0.000000 21 0.3378 0.3378
- 0.000000 22 -0.4598 -0.4598
- 0.000000 23 -0.2060 -0.2060
- 0.000000 24 0.2029 0.2029
- 0.000000 25 -0.0332 -0.0332
- 0.000000 26 -0.1518 -0.1518
- 0.000000 27 -0.1945 -0.1945
- 0.000000 28 -0.0276 -0.0276
- 0.000000 29 0.1580 0.1580
- 0.000000 30 1.5237 1.5237
- 0.000000 31 0.0115 0.0115
- 0.000000 32 -0.1825 -0.1825
- 0.000000 33 0.0115 0.0115
- 0.000000 34 1.4733 1.4733
- 0.000000 35 -0.0457 -0.0457
- 0.000000 36 -0.1825 -0.1825
- 0.000000 37 -0.0457 -0.0457
- 0.000000 38 2.4123 2.4123
- 0.050000 0 0.2922 0.2922
- 0.050000 1 0.1922 0.1922
- 0.050000 2 0.0855 0.0855
- 0.050000 3 -0.2963 -0.2963
- 0.050000 4 0.3881 0.3881
- 0.050000 5 0.0437 0.0437
- 0.050000 6 -0.3853 -0.3853
- 0.050000 7 -0.3607 -0.3607
- 0.050000 8 0.3046 0.3046
- 0.050000 9 0.4894 0.4894
- 0.050000 10 -0.4478 -0.4478
- 0.050000 11 0.0170 0.0170
- 0.050000 12 0.4238 0.4238
- 0.050000 13 0.5389 0.5389
- 0.050000 14 0.0323 0.0323
- 0.050000 15 -0.3951 -0.3951
- 0.050000 16 0.3066 0.3066
- 0.050000 17 -0.2632 -0.2632
- 0.050000 18 -0.4779 -0.4779
- 0.050000 19 -0.4842 -0.4842
- 0.050000 20 -0.0331 -0.0331
- 0.050000 21 0.3316 0.3316
- 0.050000 22 -0.4619 -0.4619
- 0.050000 23 -0.1763 -0.1763
- 0.050000 24 0.1982 0.1982
- 0.050000 25 -0.0357 -0.0357
- 0.050000 26 -0.1594 -0.1594
- 0.050000 27 -0.1805 -0.1805
- 0.050000 28 -0.0447 -0.0447
- 0.050000 29 0.1489 0.1489
- 0.050000 30 1.5486 1.5486
- 0.050000 31 0.0117 0.0117
- 0.050000 32 -0.1787 -0.1787
- 0.050000 33 0.0117 0.0117
- 0.050000 34 1.4131 1.4131
- 0.050000 35 -0.0903 -0.0903
- 0.050000 36 -0.1787 -0.1787
- 0.050000 37 -0.0903 -0.0903
- 0.050000 38 2.3847 2.3847
- 0.100000 0 0.2868 0.2868
- 0.100000 1 0.2227 0.2227
- 0.100000 2 0.1018 0.1018
- 0.100000 3 -0.3098 -0.3098
- 0.100000 4 0.3646 0.3646
- 0.100000 5 0.0698 0.0698
- 0.100000 6 -0.3693 -0.3693
- 0.100000 7 -0.3493 -0.3493
- 0.100000 8 0.3021 0.3021
- 0.100000 9 0.4580 0.4580
- 0.100000 10 -0.4736 -0.4736
- 0.100000 11 0.0569 0.0569
- 0.100000 12 0.4509 0.4509
- 0.100000 13 0.5657 0.5657
- 0.100000 14 0.0216 0.0216
- 0.100000 15 -0.3870 -0.3870
- 0.100000 16 0.3012 0.3012
- 0.100000 17 -0.3213 -0.3213
- 0.100000 18 -0.4737 -0.4737
- 0.100000 19 -0.5234 -0.5234
- 0.100000 20 -0.0355 -0.0355
- 0.100000 21 0.3428 0.3428
- 0.100000 22 -0.4281 -0.4281
- 0.100000 23 -0.2265 -0.2265
- 0.100000 24 0.1685 0.1685
- 0.100000 25 -0.1131 -0.1131
- 0.100000 26 -0.1160 -0.1160
- 0.100000 27 -0.1673 -0.1673
- 0.100000 28 -0.0498 -0.0498
- 0.100000 29 0.1471 0.1471
- 0.100000 30 1.6012 1.6012
- 0.100000 31 0.0119 0.0119
- 0.100000 32 -0.1562 -0.1562
- 0.100000 33 0.0119 0.0119
- 0.100000 34 1.3419 1.3419
- 0.100000 35 -0.0772 -0.0772
- 0.100000 36 -0.1562 -0.1562
- 0.100000 37 -0.0772 -0.0772
- 0.100000 38 2.4105 2.4105
- 0.150000 0 0.2663 0.2663
- 0.150000 1 0.2566 0.2566
- 0.150000 2 0.1199 0.1199
- 0.150000 3 -0.3179 -0.3179
- 0.150000 4 0.3869 0.3869
- 0.150000 5 0.1016 0.1016
- 0.150000 6 -0.3303 -0.3303
- 0.150000 7 -0.3223 -0.3223
- 0.150000 8 0.2965 0.2965
- 0.150000 9 0.4174 0.4174
- 0.150000 10 -0.4746 -0.4746
- 0.150000 11 0.0908 0.0908
- 0.150000 12 0.4502 0.4502
- 0.150000 13 0.5255 0.5255
- 0.150000 14 0.0136 0.0136
- 0.150000 15 -0.3639 -0.3639
- 0.150000 16 0.2792 0.2792
- 0.150000 17 -0.3581 -0.3581
- 0.150000 18 -0.4536 -0.4536
- 0.150000 19 -0.5486 -0.5486
- 0.150000 20 -0.0527 -0.0527
- 0.150000 21 0.3333 0.3333
- 0.150000 22 -0.3978 -0.3978
- 0.150000 23 -0.2670 -0.2670
- 0.150000 24 0.1374 0.1374
- 0.150000 25 -0.1500 -0.1500
- 0.150000 26 -0.0847 -0.0847
- 0.150000 27 -0.1389 -0.1389
- 0.150000 28 -0.0414 -0.0414
- 0.150000 29 0.1400 0.1400
- 0.150000 30 1.5931 1.5931
- 0.150000 31 -0.0087 -0.0087
- 0.150000 32 -0.1400 -0.1400
- 0.150000 33 -0.0087 -0.0087
- 0.150000 34 1.2403 1.2403
- 0.150000 35 -0.0378 -0.0378
- 0.150000 36 -0.1400 -0.1400
- 0.150000 37 -0.0378 -0.0378
- 0.150000 38 2.3941 2.3941
- 0.200000 0 0.2496 0.2496
- 0.200000 1 0.3129 0.3129
- 0.200000 2 0.1317 0.1317
- 0.200000 3 -0.3044 -0.3044
- 0.200000 4 0.4052 0.4052
- 0.200000 5 0.0983 0.0983
- 0.200000 6 -0.3003 -0.3003
- 0.200000 7 -0.2934 -0.2934
- 0.200000 8 0.2933 0.2933
- 0.200000 9 0.3802 0.3802
- 0.200000 10 -0.4720 -0.4720
- 0.200000 11 0.0776 0.0776
- 0.200000 12 0.4251 0.4251
- 0.200000 13 0.4943 0.4943
- 0.200000 14 0.0204 0.0204
- 0.200000 15 -0.3579 -0.3579
- 0.200000 16 0.2633 0.2633
- 0.200000 17 -0.3636 -0.3636
- 0.200000 18 -0.4105 -0.4105
- 0.200000 19 -0.5453 -0.5453
- 0.200000 20 -0.0493 -0.0493
- 0.200000 21 0.3188 0.3188
- 0.200000 22 -0.3850 -0.3850
- 0.200000 23 -0.2802 -0.2802
- 0.200000 24 0.1134 0.1134
- 0.200000 25 -0.1161 -0.1161
- 0.200000 26 -0.0641 -0.0641
- 0.200000 27 -0.1140 -0.1140
- 0.200000 28 -0.0215 -0.0215
- 0.200000 29 0.1358 0.1358
- 0.200000 30 1.5276 1.5276
- 0.200000 31 -0.0195 -0.0195
- 0.200000 32 -0.1287 -0.1287
- 0.200000 33 -0.0195 -0.0195
- 0.200000 34 1.2111 1.2111
- 0.200000 35 0.0077 0.0077
- 0.200000 36 -0.1287 -0.1287
- 0.200000 37 0.0077 0.0077
- 0.200000 38 2.2957 2.2957
+#! FIELDS time parameter s2.mean
+ 0.000000 0 0.3185
+ 0.000000 1 0.2157
+ 0.000000 2 0.1024
+ 0.000000 3 -0.3069
+ 0.000000 4 0.4357
+ 0.000000 5 0.0314
+ 0.000000 6 -0.3808
+ 0.000000 7 -0.3976
+ 0.000000 8 0.2787
+ 0.000000 9 0.5133
+ 0.000000 10 -0.4440
+ 0.000000 11 0.0141
+ 0.000000 12 0.3969
+ 0.000000 13 0.5493
+ 0.000000 14 0.0211
+ 0.000000 15 -0.3772
+ 0.000000 16 0.2829
+ 0.000000 17 -0.2163
+ 0.000000 18 -0.5101
+ 0.000000 19 -0.4685
+ 0.000000 20 -0.0317
+ 0.000000 21 0.3378
+ 0.000000 22 -0.4598
+ 0.000000 23 -0.2060
+ 0.000000 24 0.2029
+ 0.000000 25 -0.0332
+ 0.000000 26 -0.1518
+ 0.000000 27 -0.1945
+ 0.000000 28 -0.0276
+ 0.000000 29 0.1580
+ 0.000000 30 1.5237
+ 0.000000 31 0.0115
+ 0.000000 32 -0.1825
+ 0.000000 33 0.0115
+ 0.000000 34 1.4733
+ 0.000000 35 -0.0457
+ 0.000000 36 -0.1825
+ 0.000000 37 -0.0457
+ 0.000000 38 2.4123
+ 0.050000 0 0.2922
+ 0.050000 1 0.1922
+ 0.050000 2 0.0855
+ 0.050000 3 -0.2963
+ 0.050000 4 0.3881
+ 0.050000 5 0.0437
+ 0.050000 6 -0.3853
+ 0.050000 7 -0.3607
+ 0.050000 8 0.3046
+ 0.050000 9 0.4894
+ 0.050000 10 -0.4478
+ 0.050000 11 0.0170
+ 0.050000 12 0.4238
+ 0.050000 13 0.5389
+ 0.050000 14 0.0323
+ 0.050000 15 -0.3951
+ 0.050000 16 0.3066
+ 0.050000 17 -0.2632
+ 0.050000 18 -0.4779
+ 0.050000 19 -0.4842
+ 0.050000 20 -0.0331
+ 0.050000 21 0.3316
+ 0.050000 22 -0.4619
+ 0.050000 23 -0.1763
+ 0.050000 24 0.1982
+ 0.050000 25 -0.0357
+ 0.050000 26 -0.1594
+ 0.050000 27 -0.1805
+ 0.050000 28 -0.0447
+ 0.050000 29 0.1489
+ 0.050000 30 1.5486
+ 0.050000 31 0.0117
+ 0.050000 32 -0.1787
+ 0.050000 33 0.0117
+ 0.050000 34 1.4131
+ 0.050000 35 -0.0903
+ 0.050000 36 -0.1787
+ 0.050000 37 -0.0903
+ 0.050000 38 2.3847
+ 0.100000 0 0.2868
+ 0.100000 1 0.2227
+ 0.100000 2 0.1018
+ 0.100000 3 -0.3098
+ 0.100000 4 0.3646
+ 0.100000 5 0.0698
+ 0.100000 6 -0.3693
+ 0.100000 7 -0.3493
+ 0.100000 8 0.3021
+ 0.100000 9 0.4580
+ 0.100000 10 -0.4736
+ 0.100000 11 0.0569
+ 0.100000 12 0.4509
+ 0.100000 13 0.5657
+ 0.100000 14 0.0216
+ 0.100000 15 -0.3870
+ 0.100000 16 0.3012
+ 0.100000 17 -0.3213
+ 0.100000 18 -0.4737
+ 0.100000 19 -0.5234
+ 0.100000 20 -0.0355
+ 0.100000 21 0.3428
+ 0.100000 22 -0.4281
+ 0.100000 23 -0.2265
+ 0.100000 24 0.1685
+ 0.100000 25 -0.1131
+ 0.100000 26 -0.1160
+ 0.100000 27 -0.1673
+ 0.100000 28 -0.0498
+ 0.100000 29 0.1471
+ 0.100000 30 1.6012
+ 0.100000 31 0.0119
+ 0.100000 32 -0.1562
+ 0.100000 33 0.0119
+ 0.100000 34 1.3419
+ 0.100000 35 -0.0772
+ 0.100000 36 -0.1562
+ 0.100000 37 -0.0772
+ 0.100000 38 2.4105
+ 0.150000 0 0.2663
+ 0.150000 1 0.2566
+ 0.150000 2 0.1199
+ 0.150000 3 -0.3179
+ 0.150000 4 0.3869
+ 0.150000 5 0.1016
+ 0.150000 6 -0.3303
+ 0.150000 7 -0.3223
+ 0.150000 8 0.2965
+ 0.150000 9 0.4174
+ 0.150000 10 -0.4746
+ 0.150000 11 0.0908
+ 0.150000 12 0.4502
+ 0.150000 13 0.5255
+ 0.150000 14 0.0136
+ 0.150000 15 -0.3639
+ 0.150000 16 0.2792
+ 0.150000 17 -0.3581
+ 0.150000 18 -0.4536
+ 0.150000 19 -0.5486
+ 0.150000 20 -0.0527
+ 0.150000 21 0.3333
+ 0.150000 22 -0.3978
+ 0.150000 23 -0.2670
+ 0.150000 24 0.1374
+ 0.150000 25 -0.1500
+ 0.150000 26 -0.0847
+ 0.150000 27 -0.1389
+ 0.150000 28 -0.0414
+ 0.150000 29 0.1400
+ 0.150000 30 1.5931
+ 0.150000 31 -0.0087
+ 0.150000 32 -0.1400
+ 0.150000 33 -0.0087
+ 0.150000 34 1.2403
+ 0.150000 35 -0.0378
+ 0.150000 36 -0.1400
+ 0.150000 37 -0.0378
+ 0.150000 38 2.3941
+ 0.200000 0 0.2496
+ 0.200000 1 0.3129
+ 0.200000 2 0.1317
+ 0.200000 3 -0.3044
+ 0.200000 4 0.4052
+ 0.200000 5 0.0983
+ 0.200000 6 -0.3003
+ 0.200000 7 -0.2934
+ 0.200000 8 0.2933
+ 0.200000 9 0.3802
+ 0.200000 10 -0.4720
+ 0.200000 11 0.0776
+ 0.200000 12 0.4251
+ 0.200000 13 0.4943
+ 0.200000 14 0.0204
+ 0.200000 15 -0.3579
+ 0.200000 16 0.2633
+ 0.200000 17 -0.3636
+ 0.200000 18 -0.4105
+ 0.200000 19 -0.5453
+ 0.200000 20 -0.0493
+ 0.200000 21 0.3188
+ 0.200000 22 -0.3850
+ 0.200000 23 -0.2802
+ 0.200000 24 0.1134
+ 0.200000 25 -0.1161
+ 0.200000 26 -0.0641
+ 0.200000 27 -0.1140
+ 0.200000 28 -0.0215
+ 0.200000 29 0.1358
+ 0.200000 30 1.5276
+ 0.200000 31 -0.0195
+ 0.200000 32 -0.1287
+ 0.200000 33 -0.0195
+ 0.200000 34 1.2111
+ 0.200000 35 0.0077
+ 0.200000 36 -0.1287
+ 0.200000 37 0.0077
+ 0.200000 38 2.2957
diff --git a/regtest/rt22/plumed.dat b/regtest/rt22/plumed.dat
index 04735c319..2a0dd5d79 100644
--- a/regtest/rt22/plumed.dat
+++ b/regtest/rt22/plumed.dat
@@ -6,9 +6,7 @@ PRINT ARG=s1.* FILE=colvar1 FMT=%8.4f
DUMPDERIVATIVES ARG=s1.*,s1n.* FILE=derivatives1 FMT=%8.4f
COORDINATIONNUMBER SPECIES=1-10 SWITCH={RATIONAL R_0=1.0} LABEL=c1
-COORDINATIONNUMBER SPECIES=1-10 SWITCH={RATIONAL R_0=1.0} LABEL=c1n
SUBCELL ARG=c1 YLOWER=0.0 YUPPER=0.5 SIGMA=0.1 MEAN LABEL=s2
-SUBCELL ARG=c1n YLOWER=0.0 YUPPER=0.5 SIGMA=0.1 MEAN NUMERICAL_DERIVATIVES LABEL=s2n
PRINT ARG=s2.* FILE=colvar2 FMT=%8.4f
-DUMPDERIVATIVES ARG=s2.*,s2n.* FILE=derivatives2 FMT=%8.4f
+DUMPDERIVATIVES ARG=s2.* FILE=derivatives2 FMT=%8.4f
diff --git a/src/core/ActionAtomistic.cpp b/src/core/ActionAtomistic.cpp
index de9186a05..b3c860bd8 100644
--- a/src/core/ActionAtomistic.cpp
+++ b/src/core/ActionAtomistic.cpp
@@ -193,6 +193,24 @@ void ActionAtomistic::retrieveAtoms(){
if(cc && cc->checkIsEnergy()) energy=atoms.getEnergy();
}
+void ActionAtomistic::setForcesOnAtoms( const std::vector<double>& forcesToApply, unsigned ind ){
+ for(unsigned i=0;i<indexes.size();++i){
+ forces[i][0]=forcesToApply[ind]; ind++;
+ forces[i][1]=forcesToApply[ind]; ind++;
+ forces[i][2]=forcesToApply[ind]; ind++;
+ }
+ virial(0,0)=forcesToApply[ind]; ind++;
+ virial(0,1)=forcesToApply[ind]; ind++;
+ virial(0,2)=forcesToApply[ind]; ind++;
+ virial(1,0)=forcesToApply[ind]; ind++;
+ virial(1,1)=forcesToApply[ind]; ind++;
+ virial(1,2)=forcesToApply[ind]; ind++;
+ virial(2,0)=forcesToApply[ind]; ind++;
+ virial(2,1)=forcesToApply[ind]; ind++;
+ virial(2,2)=forcesToApply[ind]; ind++;
+ plumed_dbg_assert( ind==forcesToApply.size() );
+}
+
void ActionAtomistic::applyForces(){
vector<Vector> & f(atoms.forces);
Tensor & v(atoms.virial);
diff --git a/src/core/ActionAtomistic.h b/src/core/ActionAtomistic.h
index 9ea863b5f..5020fca02 100644
--- a/src/core/ActionAtomistic.h
+++ b/src/core/ActionAtomistic.h
@@ -96,6 +96,8 @@ public:
void parseAtomList(const std::string&key,const int num, std::vector<AtomNumber> &t);
/// Get reference to Pbc
const Pbc & getPbc() const;
+/// Add the forces to the atoms
+ void setForcesOnAtoms( const std::vector<double>& forcesToApply, unsigned ind=0 );
public:
// virtual functions:
diff --git a/src/multicolvar/ActionVolume.cpp b/src/multicolvar/ActionVolume.cpp
index 05e26304a..9c3de3093 100644
--- a/src/multicolvar/ActionVolume.cpp
+++ b/src/multicolvar/ActionVolume.cpp
@@ -82,8 +82,6 @@ void ActionVolume::doJobsRequiredBeforeTaskList(){
}
bool ActionVolume::performTask( const unsigned& j ){
- nderivatives=mycolv->getNumberOfDerivatives();
-
Vector catom_pos=mycolv->retrieveCentralAtomPos( derivativesOfFractionalCoordinates() );
double weight; Vector wdf;
@@ -101,32 +99,33 @@ bool ActionVolume::performTask( const unsigned& j ){
} else {
// Copy derivatives of the colvar and the value of the colvar
double colv=mycolv->getElementValue(0); setElementValue( 0, colv );
- for(unsigned i=0;i<mycolv->nderivatives;++i){
- addElementDerivative( mycolv->getOutputDerivativeIndex(mycolv->current,i), mycolv->getElementDerivative(i) );
+ for(unsigned i=mycolv->getFirstDerivativeToMerge();i<mycolv->getNumberOfDerivatives();i=mycolv->getNextDerivativeToMerge(i)){
+ addElementDerivative( i, mycolv->getElementDerivative(i) );
}
double ww=mycolv->getElementValue(1);
setElementValue( 1, ww*weight );
+ unsigned nder=getNumberOfDerivatives();
// Add derivatives of weight if we have a weight
if( mycolv->weightHasDerivatives ){
- for(unsigned i=0;i<mycolv->nderivatives;++i){
- addElementDerivative( nderivatives+mycolv->getOutputDerivativeIndex(mycolv->current,i), weight*mycolv->getElementDerivative(nderivatives+i) );
+ for(unsigned i=mycolv->getFirstDerivativeToMerge();i<mycolv->getNumberOfDerivatives();i=mycolv->getNextDerivativeToMerge(i)){
+ addElementDerivative( nder+i, weight*mycolv->getElementDerivative(nder+i) );
}
}
// Add derivatives of central atoms
for(unsigned i=0;i<mycolv->atomsWithCatomDer.getNumberActive();++i){
unsigned n=mycolv->atomsWithCatomDer[i];
- unsigned nx=nderivatives+3*linkIndex( n, mycolv->colvar_atoms[mycolv->current], mycolv->all_atoms);
+ unsigned nx=nder+3*linkIndex( n, mycolv->colvar_atoms[mycolv->current], mycolv->all_atoms);
addElementDerivative( nx+0, ww*mycolv->getCentralAtomDerivative(i, 0, wdf ) );
addElementDerivative( nx+1, ww*mycolv->getCentralAtomDerivative(i, 1, wdf ) );
addElementDerivative( nx+2, ww*mycolv->getCentralAtomDerivative(i, 2, wdf ) );
}
}
- // And compute the vessels
- return calculateAllVessels();
+ // Only continue if the weight is greater than tolerance
+ return ( weight>=getTolerance() );
}
void ActionVolume::calculateNumericalDerivatives(){
diff --git a/src/multicolvar/ActionVolume.h b/src/multicolvar/ActionVolume.h
index 34410ea80..346d416d2 100644
--- a/src/multicolvar/ActionVolume.h
+++ b/src/multicolvar/ActionVolume.h
@@ -68,8 +68,6 @@ public:
unsigned getNumberOfDerivatives(); // N.B. This is replacing the virtual function in ActionWithValue
/// Return the number of Colvars this is calculating
unsigned getNumberOfFunctionsInAction();
-/// Return the number of derivatives for a given colvar
- unsigned getNumberOfDerivatives( const unsigned& j );
/// Is the output quantity periodic
bool isPeriodic();
/// Do jobs required before tasks are undertaken
@@ -102,11 +100,6 @@ unsigned ActionVolume::getNumberOfFunctionsInAction(){
return mycolv->getNumberOfFunctionsInAction();
}
-inline
-unsigned ActionVolume::getNumberOfDerivatives( const unsigned& j ){
- return mycolv->getNumberOfDerivatives(j);
-}
-
}
}
#endif
diff --git a/src/multicolvar/MultiColvar.cpp b/src/multicolvar/MultiColvar.cpp
index b73be9095..1826bd78f 100644
--- a/src/multicolvar/MultiColvar.cpp
+++ b/src/multicolvar/MultiColvar.cpp
@@ -52,7 +52,7 @@ void MultiColvar::registerKeywords( Keywords& keys ){
keys.reserve("atoms-4","SPECIESB","this keyword is used for colvars such as the coordination number. It must appear with SPECIESA. For a full explanation see "
"the documentation for that keyword");
keys.addFlag("VERBOSE",false,"write a more detailed output");
- keys.add("optional","NL_STRIDE","the frequency with which the neighbor list should be updated. Between neighbour list update steps all quantities "
+ keys.add("hidden","NL_STRIDE","the frequency with which the neighbor list should be updated. Between neighbour list update steps all quantities "
"that contributed less than TOL at the previous neighbor list update step are ignored.");
ActionWithVessel::registerKeywords( keys );
}
@@ -90,7 +90,7 @@ void MultiColvar::addColvar( const std::vector<unsigned>& newatoms ){
newlist.addIndexToList( newatoms[i] );
}
if( verbose_output ) log.printf("\n");
- colvar_list.addIndexToList( colvar_atoms.size() );
+ taskList.addIndexToList( colvar_atoms.size() );
colvar_atoms.push_back( newlist );
}
@@ -100,14 +100,17 @@ void MultiColvar::readAtoms( int& natoms ){
if( keywords.exists("SPECIES") ) readSpeciesKeyword( natoms );
if( !readatoms ) error("No atoms have been read in");
- all_atoms.activateAll(); colvar_list.activateAll();
+ taskList.activateAll();
for(unsigned i=0;i<colvar_atoms.size();++i) colvar_atoms[i].activateAll();
- ActionAtomistic::requestAtoms( all_atoms.retrieveActiveList() );
+ // Set up dynamic array for derivatives
+ for(unsigned i=0;i<colvar_atoms[0].fullSize();++i) atoms_with_derivatives.addIndexToList( i );
// Set up stuff for central atoms
for(unsigned i=0;i<colvar_atoms[0].fullSize();++i) atomsWithCatomDer.addIndexToList( i );
atomsWithCatomDer.deactivateAll();
central_derivs.resize( colvar_atoms[0].fullSize() );
for(unsigned i=0;i<central_derivs.size();++i) central_derivs[i].zero();
+ // Resize all dynamic content
+ resizeDynamicArrays();
}
void MultiColvar::readBackboneAtoms( const std::vector<std::string>& backnames, std::vector<unsigned>& chain_lengths ){
@@ -329,27 +332,37 @@ void MultiColvar::readSpeciesKeyword( int& natoms ){
}
}
+void MultiColvar::resizeDynamicArrays(){
+ for(unsigned i=0;i<taskList.getNumberActive();++i){
+ unsigned n=taskList[i];
+ activateLinks( colvar_atoms[n], all_atoms );
+ }
+
+ all_atoms.updateActiveMembers();
+ ActionAtomistic::requestAtoms( all_atoms.retrieveActiveList() );
+ forcesToApply.resize( getNumberOfDerivatives() );
+}
+
void MultiColvar::prepare(){
updatetime=false;
if( reduceAtNextStep ){
- colvar_list.mpi_gatherActiveMembers( comm );
+ taskList.mpi_gatherActiveMembers( comm );
mpi_gatherActiveMembers( comm, colvar_atoms );
reduceAtNextStep=false; updatetime=true;
}
if( updateFreq>0 && (getStep()-lastUpdate)>=updateFreq ){
- colvar_list.activateAll();
- for(unsigned i=0;i<colvar_list.getNumberActive();++i) colvar_atoms[i].activateAll();
+ taskList.activateAll();
+ for(unsigned i=0;i<taskList.getNumberActive();++i) colvar_atoms[i].activateAll();
reduceAtNextStep=true; updatetime=true; lastUpdate=getStep();
}
if(updatetime){
- for(unsigned i=0;i<colvar_list.getNumberActive();++i) activateLinks( colvar_atoms[i], all_atoms );
- all_atoms.updateActiveMembers();
- ActionAtomistic::requestAtoms( all_atoms.retrieveActiveList() ); resizeFunctions();
+ resizeDynamicArrays();
+ resizeFunctions();
}
}
void MultiColvar::calculate(){
- runAllTasks( colvar_list.getNumberActive() );
+ runAllTasks();
}
Vector MultiColvar::getCentralAtom(){
@@ -368,15 +381,14 @@ const std::vector<Vector>& MultiColvar::getPositions(){
// Resize everything
if( pos.size()!=natoms ) pos.resize(natoms);
// Get the position
- for(unsigned i=0;i<natoms;++i) pos[i]=ActionAtomistic::getPosition( colvar_atoms[current][i] );
+ for(unsigned i=0;i<natoms;++i) pos[i]=ActionAtomistic::getPosition( getAtomIndex(i) );
posHasBeenSet=true;
}
return pos;
}
bool MultiColvar::performTask( const unsigned& j ){
- current=colvar_list[j]; // Store the number of the atom list we are currently working on
- nderivatives=3*colvar_atoms[current].getNumberActive() + 9; // Store the number of derivatives (for better performance)
+ atoms_with_derivatives.deactivateAll(); // Currently no atoms have derivatives
if( colvar_atoms[current].getNumberActive()==0 ) return false; // Do nothing if there are no active atoms in the colvar
posHasBeenSet=false; // This ensures we don't set the pos array more than once per step
@@ -388,7 +400,8 @@ bool MultiColvar::performTask( const unsigned& j ){
double vv=compute( j );
// Set the value of this element in ActionWithVessel
setElementValue( 0, vv );
-
+ // And finally gather the list of active atoms
+ atoms_with_derivatives.updateActiveMembers();
return true;
}
@@ -422,80 +435,56 @@ double MultiColvar::getCentralAtomDerivative( const unsigned& iatom, const unsig
return df[0]*central_derivs[iatom](0,jcomp) + df[1]*central_derivs[iatom](1,jcomp) + df[2]*central_derivs[iatom](2,jcomp);
}
-void MultiColvar::chainRuleForElementDerivatives( const unsigned& iout, const unsigned& ider, const unsigned& stride,
- const unsigned& off, const double& df, vesselbase::Vessel* valout ){
- plumed_dbg_assert( off<stride );
- unsigned nder=getNumberOfDerivatives();
- int thisatom, thispos, in=ider*nder, bstart=stride*(nder+1)*iout+stride+off;
- unsigned innat=colvar_atoms[current].getNumberActive();
- for(unsigned i=0;i<innat;++i){
- thisatom=linkIndex( i, colvar_atoms[current], all_atoms );
- plumed_dbg_assert( thisatom>=0 ); thispos=bstart+3*stride*thisatom;
- valout->addToBufferElement( thispos , df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( thispos+stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( thispos+2*stride, df*getElementDerivative(in) ); in++;
- }
-
- // Easy to merge the virial
- unsigned outnat=bstart+3*stride*getNumberOfAtoms();
- valout->addToBufferElement( outnat, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+2*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+3*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+4*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+5*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+6*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+7*stride, df*getElementDerivative(in) ); in++;
- valout->addToBufferElement( outnat+8*stride, df*getElementDerivative(in) ); in++;
-}
-
Vector MultiColvar::getSeparation( const Vector& vec1, const Vector& vec2 ) const {
if(usepbc){ return pbcDistance( vec1, vec2 ); }
else{ return delta( vec1, vec2 ); }
}
void MultiColvar::apply(){
- vector<Vector>& f(modifyForces());
- Tensor& v(modifyVirial());
+ getForcesFromVessels( forcesToApply );
+ setForcesOnAtoms( forcesToApply );
+}
- for(unsigned i=0;i<f.size();i++){
- f[i][0]=0.0;
- f[i][1]=0.0;
- f[i][2]=0.0;
+void MultiColvar::mergeDerivatives( const unsigned& ider, const double& df ){
+ unsigned vstart=getNumberOfDerivatives()*ider;
+ for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){
+ unsigned iatom=3*getAtomIndex( atoms_with_derivatives[i] );
+ accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++;
+ accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++;
+ accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) );
}
- v.zero();
-
- unsigned nat=getNumberOfAtoms();
- std::vector<double> forces(3*getNumberOfAtoms()+9);
-
- unsigned vstart=3*getNumberOfAtoms();
- for(int i=0;i<getNumberOfVessels();++i){
- if( (getPntrToVessel(i)->applyForce( forces )) ){
- for(unsigned j=0;j<nat;++j){
- f[j][0]+=forces[3*j+0];
- f[j][1]+=forces[3*j+1];
- f[j][2]+=forces[3*j+2];
- }
- v(0,0)+=forces[vstart+0];
- v(0,1)+=forces[vstart+1];
- v(0,2)+=forces[vstart+2];
- v(1,0)+=forces[vstart+3];
- v(1,1)+=forces[vstart+4];
- v(1,2)+=forces[vstart+5];
- v(2,0)+=forces[vstart+6];
- v(2,1)+=forces[vstart+7];
- v(2,2)+=forces[vstart+8];
- }
+ unsigned nvir=3*getNumberOfAtoms();
+ for(unsigned j=0;j<9;++j){
+ accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++;
}
}
-void MultiColvar::buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der ) const {
+void MultiColvar::clearDerivativesAfterTask( const unsigned& ider ){
+ unsigned vstart=getNumberOfDerivatives()*ider;
+ for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){
+ unsigned iatom=vstart+3*getAtomIndex( atoms_with_derivatives[i] );
+ setElementDerivative( iatom, 0.0 ); iatom++;
+ setElementDerivative( iatom, 0.0 ); iatom++;
+ setElementDerivative( iatom, 0.0 );
+ }
+ unsigned nvir=vstart+3*getNumberOfAtoms();
+ for(unsigned j=0;j<9;++j){
+ setElementDerivative( nvir, 0.0 ); nvir++;
+ }
+}
+
+void MultiColvar::buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der ){
+ // Clear old derivative indexes
+ for(unsigned i=0;i<active_der.size();++i) active_der[i].clear();
+
+ // Build indexes
active_der.resize( colvar_atoms.size() );
for(unsigned i=0;i<colvar_atoms.size();++i){
for(unsigned j=0;j<colvar_atoms[i].fullSize();++j){
- active_der[i].addIndexToList( 3*colvar_atoms[i][j] + 0 );
- active_der[i].addIndexToList( 3*colvar_atoms[i][j] + 1 );
- active_der[i].addIndexToList( 3*colvar_atoms[i][j] + 2 );
+ unsigned jatom=linkIndex( j, colvar_atoms[i], all_atoms );
+ active_der[i].addIndexToList( 3*jatom + 0 );
+ active_der[i].addIndexToList( 3*jatom + 1 );
+ active_der[i].addIndexToList( 3*jatom + 2 );
}
unsigned virs=3*getNumberOfAtoms();
for(unsigned j=0;j<9;++j){
diff --git a/src/multicolvar/MultiColvar.h b/src/multicolvar/MultiColvar.h
index 1c3b04e31..7f2e94449 100644
--- a/src/multicolvar/MultiColvar.h
+++ b/src/multicolvar/MultiColvar.h
@@ -46,6 +46,7 @@ class MultiColvar :
public vesselbase::ActionWithVessel
{
friend class ActionVolume;
+friend class MultiColvarFunction;
friend class StoreCentralAtomsVessel;
private:
bool usepbc;
@@ -60,8 +61,8 @@ private:
bool posHasBeenSet;
/// The list of all the atoms involved in the colvar
DynamicList<AtomNumber> all_atoms;
-/// A DynamicList containing the numbers of 1 - ncolvars
- DynamicList<unsigned> colvar_list;
+/// A dynamic list containing those atoms with derivatives
+ DynamicList<unsigned> atoms_with_derivatives;
/// The lists of the atoms involved in each of the individual colvars
/// note these refer to the atoms in all_atoms
std::vector< DynamicList<unsigned> > colvar_atoms;
@@ -72,12 +73,18 @@ private:
bool centralAtomDerivativesAreInFractional;
DynamicList<unsigned> atomsWithCatomDer;
std::vector<Tensor> central_derivs;
+/// Stuff used to merge derivatives
+ unsigned imerge_deriv, imerge_natoms;
+/// The forces we are going to apply to things
+ std::vector<double> forcesToApply;
/// Read in the various GROUP keywords
void readGroupsKeyword( int& natoms );
/// Read in the various SPECIES keywords
void readSpeciesKeyword( int& natoms );
/// Update the atoms request
void requestAtoms();
+/// Resize all the dynamic arrays (used at neighbor list update time and during setup)
+ void resizeDynamicArrays();
protected:
/// Are we on an update step
bool updatetime;
@@ -89,8 +96,8 @@ protected:
void readBackboneAtoms( const std::vector<std::string>& backnames, std::vector<unsigned>& chain_lengths );
/// Add a colvar to the set of colvars we are calculating (in practise just a list of atoms)
void addColvar( const std::vector<unsigned>& newatoms );
-/// Build lists of indexes of the derivatives from the colvar atoms arrays
- void buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der ) const ;
+/// Return the index of an atom
+ unsigned getAtomIndex( const unsigned& ) const ;
/// Get the separation between a pair of vectors
Vector getSeparation( const Vector& vec1, const Vector& vec2 ) const ;
/// Find out if it is time to do neighbor list update
@@ -133,16 +140,8 @@ public:
unsigned getNumberOfDerivatives(); // N.B. This is replacing the virtual function in ActionWithValue
/// Return the number of Colvars this is calculating
unsigned getNumberOfFunctionsInAction();
-/// Return the number of derivatives for a given colvar
- unsigned getNumberOfDerivatives( const unsigned& j );
/// Retrieve the position of the central atom
Vector retrieveCentralAtomPos( const bool& frac );
-/// Make sure we calculate the position of the central atom
- void useCentralAtom();
-/// Merge the derivatives
- void chainRuleForElementDerivatives( const unsigned& , const unsigned& , const unsigned& , const unsigned& , const double& , vesselbase::Vessel* );
-/// Also used for derivative merging
- unsigned getOutputDerivativeIndex( const unsigned& ival, const unsigned& i );
/// Turn of atom requests when this colvar is deactivated cos its small
void deactivate_task();
/// Turn on atom requests when the colvar is activated
@@ -153,6 +152,14 @@ public:
virtual void calculateWeight();
/// And a virtual function which actually computes the colvar
virtual double compute( const unsigned& j )=0;
+/// These are replacing the virtual methods in ActionWithVessel
+ void mergeDerivatives( const unsigned& ider, const double& df );
+ unsigned getFirstDerivativeToMerge();
+ unsigned getNextDerivativeToMerge( const unsigned& );
+/// Clear tempory data that is calculated for each task
+ void clearDerivativesAfterTask( const unsigned& );
+/// Build lists of indexes of the derivatives from the colvar atoms arrays
+ void buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der );
/// A virtual routine to get the position of the central atom - used for things like cv gradient
virtual Vector getCentralAtom();
/// Is this a density?
@@ -170,6 +177,22 @@ public:
StoreCentralAtomsVessel* getCentralAtoms();
};
+inline
+unsigned MultiColvar::getFirstDerivativeToMerge(){
+ imerge_deriv=0; imerge_natoms=atoms_with_derivatives.getNumberActive();
+ plumed_dbg_assert( colvar_atoms[current].isActive( atoms_with_derivatives[imerge_deriv] ) );
+ return 3*getAtomIndex( atoms_with_derivatives[imerge_deriv] );
+}
+
+inline
+unsigned MultiColvar::getNextDerivativeToMerge( const unsigned& j){
+ imerge_deriv++;
+ if( imerge_deriv>=3*imerge_natoms ) return 3*getNumberOfAtoms() - 3*imerge_natoms + imerge_deriv;
+ unsigned imerge_atom=std::floor( imerge_deriv / 3 );
+ plumed_dbg_assert( colvar_atoms[current].isActive( atoms_with_derivatives[imerge_atom] ) );
+ return 3*getAtomIndex( imerge_atom ) + imerge_deriv%3;
+}
+
inline
unsigned MultiColvar::getNumberOfDerivatives(){
return 3*getNumberOfAtoms()+9;
@@ -183,21 +206,10 @@ unsigned MultiColvar::getNumberOfFunctionsInAction(){
inline
void MultiColvar::deactivate_task(){
if( !reduceAtNextStep ) return; // Deactivating tasks only possible during neighbor list update
- colvar_list.deactivate(current); // Deactivate the colvar from the list
+ deactivateCurrentTask(); // Deactivate the colvar from the list
colvar_atoms[current].deactivateAll(); // Deactivate all atom requests for this colvar
}
-inline
-void MultiColvar::activateValue( const unsigned j ){
- colvar_atoms[j].activateAll();
- colvar_atoms[j].updateActiveMembers();
-}
-
-inline
-unsigned MultiColvar::getNumberOfDerivatives( const unsigned& j ){
- return 3*colvar_atoms[j].getNumberActive() + 9;
-}
-
inline
bool MultiColvar::isTimeForNeighborListUpdate() const {
return reduceAtNextStep;
@@ -206,7 +218,7 @@ bool MultiColvar::isTimeForNeighborListUpdate() const {
inline
void MultiColvar::removeAtomRequest( const unsigned& i ){
plumed_massert(reduceAtNextStep,"found removeAtomRequest but not during neighbor list step");
- colvar_atoms[current].deactivate(i);
+ colvar_atoms[current].deactivate( getAtomIndex(i) );
}
inline
@@ -220,61 +232,52 @@ unsigned MultiColvar::getNAtoms() const {
}
inline
-const Vector & MultiColvar::getPosition( unsigned iatom ) const {
+unsigned MultiColvar::getAtomIndex( const unsigned& iatom ) const {
plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- return ActionAtomistic::getPosition( colvar_atoms[current][iatom] );
+ return colvar_atoms[current][iatom];
+}
+
+inline
+const Vector & MultiColvar::getPosition( unsigned iatom ) const {
+ return ActionAtomistic::getPosition( getAtomIndex(iatom) );
}
inline
double MultiColvar::getMass(unsigned iatom ) const {
- plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- return ActionAtomistic::getMass( colvar_atoms[current][iatom] );
+ return ActionAtomistic::getMass( getAtomIndex(iatom) );
}
inline
double MultiColvar::getCharge(unsigned iatom ) const {
- plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- return ActionAtomistic::getCharge( colvar_atoms[current][iatom] );
+ return ActionAtomistic::getCharge( getAtomIndex(iatom) );
}
inline
AtomNumber MultiColvar::getAbsoluteIndex(unsigned iatom) const {
- plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- return ActionAtomistic::getAbsoluteIndex( colvar_atoms[current][iatom] );
+ return ActionAtomistic::getAbsoluteIndex( getAtomIndex(iatom) );
}
inline
void MultiColvar::addAtomsDerivatives(const int& iatom, const Vector& der){
- plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- addElementDerivative( 3*iatom+0, der[0] );
- addElementDerivative( 3*iatom+1, der[1] );
- addElementDerivative( 3*iatom+2, der[2] );
+ unsigned jatom=getAtomIndex(iatom);
+ atoms_with_derivatives.activate(iatom);
+ addElementDerivative( 3*jatom+0, der[0] );
+ addElementDerivative( 3*jatom+1, der[1] );
+ addElementDerivative( 3*jatom+2, der[2] );
}
inline
void MultiColvar::addBoxDerivatives(const Tensor& vir){
- int natoms=colvar_atoms[current].getNumberActive();
- addElementDerivative( 3*natoms+0, vir(0,0) );
- addElementDerivative( 3*natoms+1, vir(0,1) );
- addElementDerivative( 3*natoms+2, vir(0,2) );
- addElementDerivative( 3*natoms+3, vir(1,0) );
- addElementDerivative( 3*natoms+4, vir(1,1) );
- addElementDerivative( 3*natoms+5, vir(1,2) );
- addElementDerivative( 3*natoms+6, vir(2,0) );
- addElementDerivative( 3*natoms+7, vir(2,1) );
- addElementDerivative( 3*natoms+8, vir(2,2) );
-}
-
-inline
-unsigned MultiColvar::getOutputDerivativeIndex( const unsigned& ival, const unsigned& i ){
- plumed_dbg_assert( i<getNumberOfDerivatives( ival ) );
-
- if( i<3*colvar_atoms[ival].getNumberActive() ){
- unsigned inat=std::floor( i/3 );
- unsigned thisatom=linkIndex( inat, colvar_atoms[ival], all_atoms );
- return 3*thisatom + i%3;
- }
- return 3*getNumberOfAtoms() + i - 3*colvar_atoms[ival].getNumberActive();
+ unsigned nstart=3*getNumberOfAtoms();
+ addElementDerivative( nstart+0, vir(0,0) );
+ addElementDerivative( nstart+1, vir(0,1) );
+ addElementDerivative( nstart+2, vir(0,2) );
+ addElementDerivative( nstart+3, vir(1,0) );
+ addElementDerivative( nstart+4, vir(1,1) );
+ addElementDerivative( nstart+5, vir(1,2) );
+ addElementDerivative( nstart+6, vir(2,0) );
+ addElementDerivative( nstart+7, vir(2,1) );
+ addElementDerivative( nstart+8, vir(2,2) );
}
inline
@@ -289,8 +292,8 @@ void MultiColvar::setWeight( const double& weight ){
inline
void MultiColvar::addAtomsDerivativeOfWeight( const unsigned& iatom, const Vector& wder ){
- plumed_dbg_assert( iatom<colvar_atoms[current].getNumberActive() );
- int nstart = 3*getNumberOfAtoms() + 9 + 3*iatom;
+ unsigned nstart = 3*getNumberOfAtoms() + 9 + 3*getAtomIndex(iatom);
+ atoms_with_derivatives.activate(iatom);
addElementDerivative( nstart + 0, wder[0] );
addElementDerivative( nstart + 1, wder[1] );
addElementDerivative( nstart + 2, wder[2] );
@@ -298,7 +301,7 @@ void MultiColvar::addAtomsDerivativeOfWeight( const unsigned& iatom, const Vecto
inline
void MultiColvar::addBoxDerivativesOfWeight( const Tensor& vir ){
- int nstart = 3*getNumberOfAtoms() + 9 + 3*colvar_atoms[current].getNumberActive();
+ int nstart = 6*getNumberOfAtoms() + 9;
addElementDerivative( nstart+0, vir(0,0) );
addElementDerivative( nstart+1, vir(0,1) );
addElementDerivative( nstart+2, vir(0,2) );
diff --git a/src/multicolvar/MultiColvarFunction.cpp b/src/multicolvar/MultiColvarFunction.cpp
index 2d43cc08f..2ca918b7a 100644
--- a/src/multicolvar/MultiColvarFunction.cpp
+++ b/src/multicolvar/MultiColvarFunction.cpp
@@ -69,25 +69,25 @@ void MultiColvarFunction::addColvar( const std::vector<unsigned>& newatoms ){
if( colvar_atoms.size()>0 ) plumed_assert( colvar_atoms[0].fullSize()==newatoms.size() );
DynamicList<unsigned> newlist;
for(unsigned i=0;i<newatoms.size();++i) newlist.addIndexToList( newatoms[i] );
- colvar_list.addIndexToList( colvar_atoms.size() );
+ taskList.addIndexToList( colvar_atoms.size() );
colvar_atoms.push_back( newlist );
}
void MultiColvarFunction::completeSetup(){
- colvar_list.activateAll();
+ taskList.activateAll();
for(unsigned i=0;i<colvar_atoms.size();++i) colvar_atoms[i].activateAll();
}
void MultiColvarFunction::prepare(){
bool updatetime=false;
if( reduceAtNextStep ){
- colvar_list.mpi_gatherActiveMembers( comm );
+ taskList.mpi_gatherActiveMembers( comm );
mpi_gatherActiveMembers( comm, colvar_atoms );
reduceAtNextStep=false; updatetime=true;
}
if( updateFreq>0 && (getStep()-lastUpdate)>=updateFreq ){
- colvar_list.activateAll();
- for(unsigned i=0;i<colvar_list.getNumberActive();++i) colvar_atoms[i].activateAll();
+ taskList.activateAll();
+ for(unsigned i=0;i<taskList.getNumberActive();++i) colvar_atoms[i].activateAll();
reduceAtNextStep=true; updatetime=true; lastUpdate=getStep();
}
if(updatetime) resizeFunctions();
@@ -106,7 +106,7 @@ void MultiColvarFunction::calculate(){
// And recalculate
mycolv->calculate();
}
- runAllTasks( colvar_list.getNumberActive() );
+ runAllTasks();
}
void MultiColvarFunction::calculateNumericalDerivatives( ActionWithValue* a ){
@@ -114,9 +114,6 @@ void MultiColvarFunction::calculateNumericalDerivatives( ActionWithValue* a ){
}
bool MultiColvarFunction::performTask( const unsigned& j ){
- current=colvar_list[j];
- nderivatives=getNumberOfDerivatives();
-
// Compute a weight for this quantity
double weight=calculateWeight();
if( weight<getTolerance() ) return false;
@@ -128,7 +125,7 @@ bool MultiColvarFunction::performTask( const unsigned& j ){
setElementValue( 0, weight*val );
if( weightHasDerivatives ){
unsigned nder=getNumberOfDerivatives();
- for(unsigned i=0;i<nderivatives;++i){
+ for(unsigned i=mycolv->getFirstDerivativeToMerge();i<mycolv->getNumberOfDerivatives();i=mycolv->getNextDerivativeToMerge(i)){
setElementDerivative( i, weight*getElementDerivative(i) + val*getElementDerivative(nder+i) );
}
}
diff --git a/src/multicolvar/MultiColvarFunction.h b/src/multicolvar/MultiColvarFunction.h
index 125a00f7c..3259aec3c 100644
--- a/src/multicolvar/MultiColvarFunction.h
+++ b/src/multicolvar/MultiColvarFunction.h
@@ -42,8 +42,6 @@ private:
int updateFreq;
unsigned lastUpdate;
bool reduceAtNextStep;
-/// A DynamicList containing the numbers of 1 - ncolvars
- DynamicList<unsigned> colvar_list;
/// The lists of the atoms involved in each of the individual colvars
/// note these refer to the atoms in all_atoms
std::vector< DynamicList<unsigned> > colvar_atoms;
@@ -93,8 +91,6 @@ public:
void apply();
/// Get the number of derivatives for this action
unsigned getNumberOfDerivatives(); // N.B. This is replacing the virtual function in ActionWithValue
-/// Return the number of derivatives for a given colvar
- unsigned getNumberOfDerivatives( const unsigned& j );
/// Get the number of vectors we are looping over
unsigned getNumberOfFunctionsInAction();
/// Turn of atom requests when this colvar is deactivated cos its small
@@ -129,15 +125,10 @@ unsigned MultiColvarFunction::getNumberOfDerivatives(){
return mycolv->getNumberOfDerivatives();
}
-inline
-unsigned MultiColvarFunction::getNumberOfDerivatives( const unsigned& j ){
- return mycolv->getNumberOfDerivatives(j);
-}
-
inline
void MultiColvarFunction::deactivate_task(){
if( !reduceAtNextStep ) return; // Deactivating tasks only possible during neighbor list update
- colvar_list.deactivate(current); // Deactivate the colvar from the list
+ deactivateCurrentTask(); // Deactivate the colvar from the list
colvar_atoms[current].deactivateAll(); // Deactivate all atom requests for this colvar
}
diff --git a/src/multicolvar/StoreCentralAtomsVessel.cpp b/src/multicolvar/StoreCentralAtomsVessel.cpp
index 17c593954..3608d3c72 100644
--- a/src/multicolvar/StoreCentralAtomsVessel.cpp
+++ b/src/multicolvar/StoreCentralAtomsVessel.cpp
@@ -36,8 +36,6 @@ nspace(1 + 3*CATOMS_MAXATOMS)
{
mycolv=dynamic_cast<MultiColvar*>( getAction() );
plumed_assert( mycolv );
-
- mycolv->buildDerivativeIndexArrays( active_der );
}
void StoreCentralAtomsVessel::resize(){
@@ -50,6 +48,8 @@ void StoreCentralAtomsVessel::resize(){
start[nfunc]=bsize;
resizeBuffer( bsize );
forces.resize( mycolv->getNumberOfDerivatives() );
+ // Update the active_derivative lists
+ mycolv->buildDerivativeIndexArrays( active_der );
}
void StoreCentralAtomsVessel::prepare(){
diff --git a/src/vesselbase/ActionWithVessel.cpp b/src/vesselbase/ActionWithVessel.cpp
index 3004b3df8..67890be96 100644
--- a/src/vesselbase/ActionWithVessel.cpp
+++ b/src/vesselbase/ActionWithVessel.cpp
@@ -70,7 +70,7 @@ void ActionWithVessel::addBridgingVessel( ActionWithVessel* tome, BridgeVessel*
bv=new BridgeVessel(da);
bv->setOutputAction( tome );
functions.push_back( dynamic_cast<Vessel*>(bv) );
- resizeFunctions();
+ resizeFunctions();
}
void ActionWithVessel::readVesselKeywords(){
@@ -119,10 +119,17 @@ void ActionWithVessel::resizeFunctions(){
plumed_massert( tmpnval>1 , "There should always be at least two terms - one for the value and one for the weight");
if(tmpnval>nvals) nvals=tmpnval;
}
- nderivatives=getNumberOfDerivatives();
+ unsigned nder=getNumberOfDerivatives();
thisval.resize( nvals ); thisval_wasset.resize( nvals, false );
- derivatives.resize( nvals*nderivatives, 0.0 );
+ derivatives.resize( nvals*nder, 0.0 );
buffer.resize( bufsize );
+ tmpforces.resize( getNumberOfDerivatives() );
+}
+
+void ActionWithVessel::deactivateCurrentTask(){
+ for(unsigned i=0;i<taskList.fullSize();++i){
+ if( taskList(i)==current ) taskList.deactivate(i);
+ }
}
//Vessel* ActionWithVessel::getVessel( const std::string& name ){
@@ -139,10 +146,11 @@ void ActionWithVessel::resizeFunctions(){
void ActionWithVessel::doJobsRequiredBeforeTaskList(){
// Clear all data from previous calculations
buffer.assign(buffer.size(),0.0);
+ // Do any preparatory stuff for functions
for(unsigned j=0;j<functions.size();++j) functions[j]->prepare();
}
-void ActionWithVessel::runAllTasks( const unsigned& ntasks ){
+void ActionWithVessel::runAllTasks(){
plumed_massert( read, "you must have a call to readVesselKeywords somewhere" );
unsigned stride=comm.Get_size();
unsigned rank=comm.Get_rank();
@@ -152,7 +160,9 @@ void ActionWithVessel::runAllTasks( const unsigned& ntasks ){
// Make sure jobs are done
doJobsRequiredBeforeTaskList();
- for(unsigned i=rank;i<ntasks;i+=stride){
+ for(unsigned i=rank;i<taskList.getNumberActive();i+=stride){
+ // Store the task we are currently working on
+ current=taskList[i];
// Calculate the stuff in the loop for this action
bool keep=performTask(i);
@@ -160,9 +170,8 @@ void ActionWithVessel::runAllTasks( const unsigned& ntasks ){
// the condition is that the weight of the contribution is low
if( !keep ){
plumed_dbg_assert( thisval[1]<getTolerance() && !thisval_wasset[0] );
- // Clear everything from weight ready for next loop
- thisval_wasset[1]=false;
- for(unsigned j=0;j<nderivatives;++j) derivatives[nder+j]=0.0;
+ // Clear the derivatives
+ clearAfterTask();
// Deactivate task
deactivate_task();
continue;
@@ -176,6 +185,19 @@ void ActionWithVessel::runAllTasks( const unsigned& ntasks ){
finishComputations();
}
+void ActionWithVessel::clearAfterTask(){
+ // Clear the derivatives from this step
+ for(unsigned k=0;k<thisval.size();++k){
+ thisval_wasset[k]=false;
+ clearDerivativesAfterTask(k);
+ }
+}
+
+void ActionWithVessel::clearDerivativesAfterTask( const unsigned& ider ){
+ unsigned kstart=ider*getNumberOfDerivatives();
+ for(unsigned j=getFirstDerivativeToMerge();j<getNumberOfDerivatives();j=getNextDerivativeToMerge(j)) derivatives[ kstart+j ]=0.0;
+}
+
bool ActionWithVessel::calculateAllVessels(){
bool keep=false;
for(unsigned j=0;j<functions.size();++j){
@@ -183,12 +205,7 @@ bool ActionWithVessel::calculateAllVessels(){
// quantity is contributing more than the tolerance
if( functions[j]->calculate() ) keep=true;
}
- // Clear the derivatives from this step
- for(unsigned k=0;k<thisval.size();++k){
- thisval_wasset[k]=false;
- unsigned kstart=k*getNumberOfDerivatives();
- for(unsigned j=0;j<nderivatives;++j) derivatives[kstart+j]=0.0;
- }
+ clearAfterTask();
return keep;
}
@@ -201,14 +218,45 @@ void ActionWithVessel::finishComputations(){
}
void ActionWithVessel::chainRuleForElementDerivatives( const unsigned& iout, const unsigned& ider, const double& df, Vessel* valout ){
- chainRuleForElementDerivatives(iout,ider,1,0,df,valout);
+ current_buffer_stride=1;
+ current_buffer_start=valout->bufstart + (getNumberOfDerivatives()+1)*iout + 1;
+ mergeDerivatives( ider, df );
}
void ActionWithVessel::chainRuleForElementDerivatives( const unsigned& iout, const unsigned& ider, const unsigned& stride,
const unsigned& off, const double& df, Vessel* valout ){
plumed_dbg_assert( off<stride );
- unsigned nder=getNumberOfDerivatives(), bstart=stride*(nder+1)*iout+stride+off, vstart=nder*ider;
- for(unsigned i=0;i<nder;++i) valout->addToBufferElement( bstart+i*stride, df*derivatives[vstart+i] );
+ current_buffer_stride=stride;
+ current_buffer_start=valout->bufstart + stride*(getNumberOfDerivatives()+1)*iout + stride + off;
+ mergeDerivatives( ider, df );
+}
+
+void ActionWithVessel::mergeDerivatives( const unsigned& ider, const double& df ){
+ unsigned nder=getNumberOfDerivatives(), vstart=nder*ider;
+ for(unsigned i=getFirstDerivativeToMerge();i<getNumberOfDerivatives();i=getNextDerivativeToMerge(i)){
+ accumulateDerivative( i, df*derivatives[vstart+i] );
+ }
+}
+
+void ActionWithVessel::buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der ){
+ // Clear old derivative indexes
+ for(unsigned i=0;i<active_der.size();++i) active_der[i].clear();
+
+ // Build indexes
+ active_der.resize( getNumberOfFunctionsInAction() );
+ for(unsigned i=0;i<active_der.size();++i){
+ for(unsigned j=0;j<getNumberOfDerivatives();++j) active_der[i].addIndexToList( j );
+ }
+}
+
+void ActionWithVessel::getForcesFromVessels( std::vector<double>& forcesToApply ){
+ plumed_dbg_assert( forcesToApply.size()==getNumberOfDerivatives() );
+ forcesToApply.assign( forcesToApply.size(),0.0 );
+ for(int i=0;i<getNumberOfVessels();++i){
+ if( (functions[i]->applyForce( tmpforces )) ){
+ for(unsigned j=0;j<forcesToApply.size();++j) forcesToApply[j]+=tmpforces[j];
+ }
+ }
}
void ActionWithVessel::retrieveDomain( std::string& min, std::string& max ){
diff --git a/src/vesselbase/ActionWithVessel.h b/src/vesselbase/ActionWithVessel.h
index 3edc570c8..662a21d91 100644
--- a/src/vesselbase/ActionWithVessel.h
+++ b/src/vesselbase/ActionWithVessel.h
@@ -25,6 +25,7 @@
#include "core/ActionWithValue.h"
#include "core/ActionAtomistic.h"
#include "tools/Exception.h"
+#include "tools/DynamicList.h"
#include <vector>
namespace PLMD{
@@ -61,9 +62,13 @@ private:
/// Vector of derivatives for the object
std::vector<double> derivatives;
/// The buffers we use for mpi summing DistributionFunction objects
+ unsigned current_buffer_start;
+ unsigned current_buffer_stride;
std::vector<double> buffer;
/// Pointers to the functions we are using on each value
std::vector<Vessel*> functions;
+/// Tempory storage for forces
+ std::vector<double> tmpforces;
/// Avoid hiding of base class virtual function
using Action::deactivate;
protected:
@@ -71,8 +76,8 @@ protected:
bool weightHasDerivatives;
/// The numerical index of the task we are curently performing
unsigned current;
-/// The number of derivatives involved in the current task
- unsigned nderivatives;
+/// The list of tasks we have to perform
+ DynamicList<unsigned> taskList;
/// Add a vessel to the list of vessels
void addVessel( const std::string& name, const std::string& input, const int numlab=0, const std::string thislab="" );
void addVessel( Vessel* vv );
@@ -87,7 +92,9 @@ protected:
/// Get a pointer to the ith vessel
Vessel* getPntrToVessel( const unsigned& i );
/// Calculate the values of all the vessels
- void runAllTasks( const unsigned& ntasks );
+ void runAllTasks();
+/// Deactivate the task we are currently working on
+ void deactivateCurrentTask();
/// Finish running all the calculations
void finishComputations();
/// Resize all the functions when the number of derivatives change
@@ -97,6 +104,12 @@ protected:
/// This loops over all the vessels calculating them and also
/// sets all the element derivatives equal to zero
bool calculateAllVessels();
+/// Retrieve the forces from all the vessels (used in apply)
+ void getForcesFromVessels( std::vector<double>& forcesToApply );
+/// This is used to accumulate the derivatives when we merge using chainRuleForElementDerivatives
+ void accumulateDerivative( const unsigned& ider, const double& df );
+/// Clear tempory data that is calculated for each task
+ void clearAfterTask();
public:
static void registerKeywords(Keywords& keys);
ActionWithVessel(const ActionOptions&ao);
@@ -106,8 +119,12 @@ public:
virtual void deactivate_task()=0;
/// Merge the derivatives
void chainRuleForElementDerivatives( const unsigned&, const unsigned&, const double& , Vessel* );
- virtual void chainRuleForElementDerivatives( const unsigned&, const unsigned& , const unsigned& , const unsigned& , const double& , Vessel* );
- virtual unsigned getOutputDerivativeIndex( const unsigned& ival, const unsigned& i ){ return i; }
+ void chainRuleForElementDerivatives( const unsigned&, const unsigned& , const unsigned& , const unsigned& , const double& , Vessel* );
+ virtual void mergeDerivatives( const unsigned& ider, const double& df );
+ virtual unsigned getFirstDerivativeToMerge();
+ virtual unsigned getNextDerivativeToMerge( const unsigned& );
+ virtual void buildDerivativeIndexArrays( std::vector< DynamicList<unsigned> >& active_der );
+ virtual void clearDerivativesAfterTask( const unsigned& );
/// Are the base quantities periodic
virtual bool isPeriodic()=0;
/// What are the domains of the base quantities
@@ -116,8 +133,6 @@ public:
virtual unsigned getNumberOfFunctionsInAction()=0;
/// Get the number of derivatives for final calculated quantity
virtual unsigned getNumberOfDerivatives()=0;
-/// Get number of derivatives for ith function
- virtual unsigned getNumberOfDerivatives( const unsigned& i );
/// Do any jobs that are required before the task list is undertaken
virtual void doJobsRequiredBeforeTaskList();
/// Calculate one of the functions in the distribution
@@ -165,8 +180,13 @@ void ActionWithVessel::setElementValue( const unsigned& ival, const double& val
}
inline
-unsigned ActionWithVessel::getNumberOfDerivatives( const unsigned& i ){
- return derivatives.size();
+unsigned ActionWithVessel::getFirstDerivativeToMerge(){
+ return 0;
+}
+
+inline
+unsigned ActionWithVessel::getNextDerivativeToMerge( const unsigned& ider){
+ return ider+1;
}
inline
@@ -195,6 +215,12 @@ void ActionWithVessel::setElementDerivative( const unsigned& ider, const double&
derivatives[ider] = der;
}
+inline
+void ActionWithVessel::accumulateDerivative( const unsigned& ider, const double& der ){
+ plumed_dbg_assert( ider<getNumberOfDerivatives() );
+ buffer[current_buffer_start + current_buffer_stride*ider] += der;
+}
+
}
}
diff --git a/src/vesselbase/BridgeVessel.cpp b/src/vesselbase/BridgeVessel.cpp
index c973f6b5d..bb0715c03 100644
--- a/src/vesselbase/BridgeVessel.cpp
+++ b/src/vesselbase/BridgeVessel.cpp
@@ -55,7 +55,12 @@ void BridgeVessel::prepare(){
bool BridgeVessel::calculate(){
unsigned j;
- return myOutputAction->performTask(j);
+ bool res=myOutputAction->performTask(j);
+ if( !res ){
+ myOutputAction->clearAfterTask();
+ return false;
+ }
+ return myOutputAction->calculateAllVessels();
}
void BridgeVessel::finish(){
diff --git a/src/vesselbase/Moments.cpp b/src/vesselbase/Moments.cpp
index 924cdacf8..1f9a04cb2 100644
--- a/src/vesselbase/Moments.cpp
+++ b/src/vesselbase/Moments.cpp
@@ -39,7 +39,7 @@ public:
static void reserveKeyword( Keywords& keys );
Moments( const VesselOptions& da );
std::string description();
- void finish();
+ void performCalculationUsingAllValues();
void local_resizing();
bool applyForce( std::vector<double>& forces );
};
@@ -92,7 +92,7 @@ std::string Moments::description(){
return descri;
}
-void Moments::finish(){
+void Moments::performCalculationUsingAllValues(){
const double pi=3.141592653589793238462643383279502884197169399375105820974944592307;
unsigned nvals=getAction()->getNumberOfFunctionsInAction();
diff --git a/src/vesselbase/StoreValuesVessel.cpp b/src/vesselbase/StoreValuesVessel.cpp
index 35b79808a..50e33dd14 100644
--- a/src/vesselbase/StoreValuesVessel.cpp
+++ b/src/vesselbase/StoreValuesVessel.cpp
@@ -22,6 +22,8 @@
#include "ActionWithVessel.h"
#include "StoreValuesVessel.h"
+#define MAXDERIVATIVES 300
+
namespace PLMD {
namespace vesselbase{
@@ -42,39 +44,50 @@ void StoreValuesVessel::resize(){
bufsize=0; start.resize( nfunc +1 );
for(unsigned i=0;i<nfunc;++i){
start[i] = bufsize;
- bufsize += 1 + aa->getNumberOfDerivatives(i);
+ bufsize += 1 + MAXDERIVATIVES;
}
start[nfunc]=bufsize;
resizeBuffer( 2*bufsize ); local_resizing();
+ // Build the arrays of indexes
+ getAction()->buildDerivativeIndexArrays( active_der );
}
bool StoreValuesVessel::calculate(){
ActionWithVessel* aa=getAction();
unsigned ibuf=start[aa->current]; setBufferElement( ibuf, aa->getElementValue(0) ); ibuf++;
- for(unsigned ider=0;ider<getAction()->nderivatives;++ider){ setBufferElement( ibuf, aa->getElementDerivative(ider) ); ibuf++; }
- plumed_dbg_assert( ibuf==start[aa->current+1] );
+ for(unsigned ider=getAction()->getFirstDerivativeToMerge();ider<getAction()->getNumberOfDerivatives();ider=getAction()->getNextDerivativeToMerge(ider)){
+ active_der[aa->current].activate(ider);
+ setBufferElement( ibuf, aa->getElementDerivative(ider) ); ibuf++;
+ }
ibuf=bufsize+start[aa->current]; setBufferElement( ibuf, aa->getElementValue(1) ); ibuf++;
if(diffweight){
unsigned nder=aa->getNumberOfDerivatives();
- for(unsigned ider=0;ider<getAction()->nderivatives;++ider){ setBufferElement( ibuf, aa->getElementDerivative(nder+ider) ); ibuf++; }
- plumed_dbg_assert( ibuf==(bufsize+start[aa->current+1]) );
+ for(unsigned ider=getAction()->getFirstDerivativeToMerge();ider<getAction()->getNumberOfDerivatives();ider=getAction()->getNextDerivativeToMerge(ider)){
+ active_der[aa->current].activate(ider);
+ setBufferElement( ibuf, aa->getElementDerivative(nder+ider) ); ibuf++;
+ }
}
return true;
}
+void StoreValuesVessel::finish(){
+ mpi_gatherActiveMembers( comm, active_der );
+ performCalculationUsingAllValues();
+}
+
void StoreValuesVessel::addDerivatives( const unsigned& ival, double& pref, Value* value_out ){
- unsigned j=0; ActionWithVessel* aa=getAction();
- for(unsigned i=start[ival]+1;i<start[ival+1];++i){
- value_out->addDerivative( aa->getOutputDerivativeIndex( ival, j ), pref*getBufferElement(i) ); j++;
+ unsigned jbuf=start[ival]+1;
+ for(unsigned i=0;i<active_der[ival].getNumberActive();++i){
+ value_out->addDerivative( active_der[ival][i], pref*getBufferElement(jbuf) ); jbuf++;
}
}
void StoreValuesVessel::addWeightDerivatives( const unsigned& ival, double& pref, Value* value_out ){
if(!diffweight) return;
- unsigned j=0; ActionWithVessel* aa=getAction();
- for(unsigned i=bufsize+start[ival]+1;i<bufsize+start[ival+1];++i){
- value_out->addDerivative( aa->getOutputDerivativeIndex( ival, j ), pref*getBufferElement(i) ); j++;
+ unsigned jbuf=bufsize+start[ival]+1;
+ for(unsigned i=0;i<active_der[ival].getNumberActive();++i){
+ value_out->addDerivative( active_der[ival][i], pref*getBufferElement(jbuf) ); jbuf++;
}
}
diff --git a/src/vesselbase/StoreValuesVessel.h b/src/vesselbase/StoreValuesVessel.h
index d474d13e6..64ae4c342 100644
--- a/src/vesselbase/StoreValuesVessel.h
+++ b/src/vesselbase/StoreValuesVessel.h
@@ -26,6 +26,7 @@
#include <cstring>
#include <vector>
#include "Vessel.h"
+#include "tools/DynamicList.h"
namespace PLMD {
namespace vesselbase{
@@ -34,6 +35,7 @@ class StoreValuesVessel : public Vessel {
private:
unsigned bufsize;
std::vector<unsigned> start;
+ std::vector< DynamicList<unsigned> > active_der;
protected:
/// Are the weights differentiable
bool diffweight;
@@ -53,10 +55,14 @@ public:
unsigned getNumberOfTerms(){ return 2; }
/// This does the resizing of the buffer
void resize();
-/// This makes sure all values are stored
- bool calculate();
/// This makes sure things further down the chain are resized
virtual void local_resizing()=0;
+/// This makes sure all values are stored
+ bool calculate();
+/// Makes sure the dynamic lists are gathered
+ void finish();
+/// Do the calculation
+ virtual void performCalculationUsingAllValues()=0;
};
inline
diff --git a/src/vesselbase/Vessel.cpp b/src/vesselbase/Vessel.cpp
index 085464154..1f0a31da1 100644
--- a/src/vesselbase/Vessel.cpp
+++ b/src/vesselbase/Vessel.cpp
@@ -59,8 +59,8 @@ myname(da.myname),
numlab(da.numlab),
action(da.action),
keywords(da.keywords),
-comm(da.action->comm),
finished_read(false),
+comm(da.action->comm),
log((da.action)->log)
{
line=Tools::getWords( da.parameters );
--
GitLab