From 8080dcb547210f5d6eb672e80f7ccd731c488d57 Mon Sep 17 00:00:00 2001 From: Gareth Aneurin Tribello Date: Fri, 6 Sep 2024 17:42:58 +0100 Subject: [PATCH] Fixed to ensure that regtests gridtools/rt-spherical-integral-2, symfunc/rt-two-filter and crystdistrib/rt-bops all work with new force scheme --- .../gridtools/rt-spherical-integral-2/config | 14 ++++++------- regtest/symfunc/rt-two-filter/config | 8 ++++++++ regtest/symfunc/rt-two-filter/plumed.dat | 2 +- .../QuaternionBondProductMatrix.cpp | 4 ++-- src/function/FunctionOfMatrix.h | 3 +-- src/matrixtools/OuterProduct.cpp | 20 ++++++++++++++++--- 6 files changed, 36 insertions(+), 15 deletions(-) diff --git a/regtest/gridtools/rt-spherical-integral-2/config b/regtest/gridtools/rt-spherical-integral-2/config index 704f414c67..970da428cf 100644 --- a/regtest/gridtools/rt-spherical-integral-2/config +++ b/regtest/gridtools/rt-spherical-integral-2/config @@ -3,10 +3,10 @@ plumed_modules=adjmat # this is to test a different name arg="--plumed plumed.dat --ixyz trajectory.xyz --dump-forces forces --dump-forces-fmt %8.4f" # --debug-forces forces.num" -# function plumed_regtest_before() { -# export PLUMED_FORBID_CHAINS=yes -# } -# -# function plumed_regtest_after() { -# export PLUMED_FORBID_CHAINS=no -# } +function plumed_regtest_before() { + export PLUMED_FORBID_CHAINS=yes +} + +function plumed_regtest_after() { + export PLUMED_FORBID_CHAINS=no +} diff --git a/regtest/symfunc/rt-two-filter/config b/regtest/symfunc/rt-two-filter/config index 8d536cf0b8..3550337ee3 100644 --- a/regtest/symfunc/rt-two-filter/config +++ b/regtest/symfunc/rt-two-filter/config @@ -2,3 +2,11 @@ type=driver plumed_modules="symfunc adjmat" # this is to test a different name arg="--plumed plumed.dat --ixyz gas-one.xyz --dump-forces forces --dump-forces-fmt=%8.4f" #--debug-forces forces.num" + +function plumed_regtest_before() { + export PLUMED_FORBID_CHAINS=yes +} + +function plumed_regtest_after() { + export PLUMED_FORBID_CHAINS=no +} diff --git a/regtest/symfunc/rt-two-filter/plumed.dat b/regtest/symfunc/rt-two-filter/plumed.dat index 8d2c31f8fc..a8a5dd3f8f 100644 --- a/regtest/symfunc/rt-two-filter/plumed.dat +++ b/regtest/symfunc/rt-two-filter/plumed.dat @@ -2,7 +2,7 @@ c1: COORDINATIONNUMBER SPECIES=1-150 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0} c1m: MORE_THAN ARG=c1 SWITCH={RATIONAL D_0=2.0 R_0=0.1} con: CONTACT_MATRIX GROUP=c1 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0} -dot: OUTER_PRODUCT ARG=c1m,c1m +dot: OUTER_PRODUCT ARG=c1m,c1m MASK=con mm: MATHEVAL ARG=con,dot FUNC=x*y PERIODIC=NO ones: ONES SIZE=150 diff --git a/src/crystdistrib/QuaternionBondProductMatrix.cpp b/src/crystdistrib/QuaternionBondProductMatrix.cpp index 763180f966..e3e77ed866 100644 --- a/src/crystdistrib/QuaternionBondProductMatrix.cpp +++ b/src/crystdistrib/QuaternionBondProductMatrix.cpp @@ -127,9 +127,9 @@ unsigned QuaternionBondProductMatrix::getNumberOfColumns() const { } void QuaternionBondProductMatrix::setupForTask( const unsigned& task_index, std::vector& indices, MultiValue& myvals ) const { - unsigned start_n = getPntrToArgument(4)->getShape()[0], size_v = getPntrToArgument(4)->getShape()[1]; + Value* myarg = getPntrToArgument(4); unsigned start_n = myarg->getShape()[0], size_v = myarg->getRowLength(task_index); if( indices.size()!=size_v+1 ) indices.resize( size_v+1 ); - for(unsigned i=0; igetRowIndex( task_index, i); myvals.setSplitIndex( size_v + 1 ); } diff --git a/src/function/FunctionOfMatrix.h b/src/function/FunctionOfMatrix.h index fe866b5d16..1e7d32a013 100644 --- a/src/function/FunctionOfMatrix.h +++ b/src/function/FunctionOfMatrix.h @@ -343,8 +343,7 @@ void FunctionOfMatrix::performTask( const std::string& controller, const unsi if( getPntrToArgument(j)->getRank()==2 ) { for(int i=0; igetPositionInStream(); - unsigned myind = base + getPntrToArgument(j)->getNumberOfColumns()*index1 + ind2; - if( getConstPntrToComponent(i)->valueIsStored() && getConstPntrToComponent(i)->getNumberOfColumns()getShape()[1] ) myind = base + getPntrToArgument(j)->getNumberOfColumns()*index1 + getConstPntrToComponent(i)->getRowLength(index1); + unsigned myind = base + getPntrToArgument(j)->getIndexInStore( index1*getPntrToArgument(j)->getShape()[1] + ind2 ); myvals.addDerivative( ostrn, myind, derivatives(i,j) ); myvals.updateIndex( ostrn, myind ); } diff --git a/src/matrixtools/OuterProduct.cpp b/src/matrixtools/OuterProduct.cpp index 2a15a362ce..d5a46381c6 100644 --- a/src/matrixtools/OuterProduct.cpp +++ b/src/matrixtools/OuterProduct.cpp @@ -55,7 +55,7 @@ class OuterProduct : public ActionWithMatrix { PLUMED_REGISTER_ACTION(OuterProduct,"OUTER_PRODUCT") void OuterProduct::registerKeywords( Keywords& keys ) { - ActionWithMatrix::registerKeywords(keys); keys.use("ARG"); + ActionWithMatrix::registerKeywords(keys); keys.use("ARG"); keys.use("MASK"); keys.add("compulsory","FUNC","x*y","the function of the input vectors that should be put in the elements of the outer product"); keys.addFlag("ELEMENTS_ON_DIAGONAL_ARE_ZERO",false,"set all diagonal elements to zero"); keys.setValueDescription("a matrix containing the outer product of the two input vectors that was obtained using the function that was input"); @@ -67,9 +67,15 @@ OuterProduct::OuterProduct(const ActionOptions&ao): domin(false), domax(false) { - if( getNumberOfArguments()!=2 ) error("should be two arguments to this action, a matrix and a vector"); + unsigned nargs=getNumberOfArguments(); if( getNumberOfMasks()>0 ) nargs = nargs - getNumberOfMasks(); + if( nargs!=2 ) error("should be two arguments to this action, they should both be vectors"); if( getPntrToArgument(0)->getRank()!=1 || getPntrToArgument(0)->hasDerivatives() ) error("first argument to this action should be a vector"); if( getPntrToArgument(1)->getRank()!=1 || getPntrToArgument(1)->hasDerivatives() ) error("first argument to this action should be a vector"); + if( getNumberOfMasks()==1 ) { + if( getPntrToArgument(2)->getRank()!=2 || getPntrToArgument(2)->hasDerivatives() ) error("mask argument should be a matrix"); + if( getPntrToArgument(2)->getShape()[0]!=getPntrToArgument(0)->getShape()[0] ) error("mask argument has wrong size"); + if( getPntrToArgument(2)->getShape()[1]!=getPntrToArgument(1)->getShape()[0] ) error("mask argument has wrong size"); + } std::string func; parse("FUNC",func); if( func=="min") { @@ -106,7 +112,15 @@ void OuterProduct::prepare() { } void OuterProduct::setupForTask( const unsigned& task_index, std::vector& indices, MultiValue& myvals ) const { - unsigned start_n = getPntrToArgument(0)->getShape()[0], size_v = getPntrToArgument(1)->getShape()[0]; + unsigned start_n = getPntrToArgument(0)->getShape()[0]; + if( getNumberOfMasks()>0 ) { + Value* maskarg = getPntrToArgument(2); unsigned size_v = maskarg->getRowLength(task_index); + if( indices.size()!=size_v+1 ) indices.resize( size_v+1 ); + for(unsigned i=0; igetRowIndex( task_index, i ); + myvals.setSplitIndex( 1 + size_v ); return; + } + + unsigned size_v = getPntrToArgument(1)->getShape()[0]; if( diagzero ) { if( indices.size()!=size_v ) indices.resize( size_v ); unsigned k=1;