Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 : Copyright (c) 2013-2023 The plumed team 3 : (see the PEOPLE file at the root of the distribution for a list of names) 4 : 5 : See http://www.plumed.org for more information. 6 : 7 : This file is part of plumed, version 2. 8 : 9 : plumed is free software: you can redistribute it and/or modify 10 : it under the terms of the GNU Lesser General Public License as published by 11 : the Free Software Foundation, either version 3 of the License, or 12 : (at your option) any later version. 13 : 14 : plumed is distributed in the hope that it will be useful, 15 : but WITHOUT ANY WARRANTY; without even the implied warranty of 16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 17 : GNU Lesser General Public License for more details. 18 : 19 : You should have received a copy of the GNU Lesser General Public License 20 : along with plumed. If not, see <http://www.gnu.org/licenses/>. 21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */ 22 : #include "RMSDBase.h" 23 : #include "MetricRegister.h" 24 : #include "tools/RMSD.h" 25 : 26 : namespace PLMD { 27 : 28 : class SimpleRMSD : public RMSDBase { 29 : private: 30 : RMSD myrmsd; 31 : public: 32 : explicit SimpleRMSD( const ReferenceConfigurationOptions& ro ); 33 : void read( const PDB& ) override; 34 : double calc( const std::vector<Vector>& pos, ReferenceValuePack& myder, const bool& squared ) const override; 35 4 : bool pcaIsEnabledForThisReference() override { return true; } 36 68 : void setupPCAStorage( ReferenceValuePack& mypack ) override { 37 68 : mypack.switchOnPCAOption(); mypack.getAtomsDisplacementVector().resize( getNumberOfAtoms() ); 38 68 : } 39 : void extractAtomicDisplacement( const std::vector<Vector>& pos, std::vector<Vector>& direction ) const override; 40 : double projectAtomicDisplacementOnVector( const bool& normalized, const std::vector<Vector>& vecs, ReferenceValuePack& mypack ) const override; 41 : }; 42 : 43 10477 : PLUMED_REGISTER_METRIC(SimpleRMSD,"SIMPLE") 44 : 45 29 : SimpleRMSD::SimpleRMSD( const ReferenceConfigurationOptions& ro ): 46 : ReferenceConfiguration( ro ), 47 29 : RMSDBase( ro ) 48 : { 49 29 : } 50 : 51 25 : void SimpleRMSD::read( const PDB& pdb ) { 52 25 : readReference( pdb ); 53 25 : } 54 : 55 160 : double SimpleRMSD::calc( const std::vector<Vector>& pos, ReferenceValuePack& myder, const bool& squared ) const { 56 160 : if( myder.getAtomsDisplacementVector().size()!=pos.size() ) myder.getAtomsDisplacementVector().resize( pos.size() ); 57 160 : double d=myrmsd.simpleAlignment( getAlign(), getDisplace(), pos, getReferencePositions(), myder.getAtomVector(), myder.getAtomsDisplacementVector(), squared ); 58 2509 : myder.clear(); for(unsigned i=0; i<pos.size(); ++i) myder.setAtomDerivatives( i, myder.getAtomVector()[i] ); 59 160 : if( !myder.updateComplete() ) myder.updateDynamicLists(); 60 160 : return d; 61 : } 62 : 63 0 : void SimpleRMSD::extractAtomicDisplacement( const std::vector<Vector>& pos, std::vector<Vector>& direction ) const { 64 0 : std::vector<Vector> tder( getNumberOfAtoms() ); 65 0 : myrmsd.simpleAlignment( getAlign(), getDisplace(), pos, getReferencePositions(), tder, direction, true ); 66 0 : for(unsigned i=0; i<pos.size(); ++i) direction[i] = getDisplace()[i]*direction[i]; 67 0 : } 68 : 69 88 : double SimpleRMSD::projectAtomicDisplacementOnVector( const bool& normalized, const std::vector<Vector>& vecs, ReferenceValuePack& mypack ) const { 70 88 : plumed_dbg_assert( mypack.calcUsingPCAOption() ); Vector comder; comder.zero(); 71 550 : for(unsigned j=0; j<vecs.size(); ++j) { 72 1848 : for(unsigned k=0; k<3; ++k) comder[k] += getAlign()[j]*vecs[j][k]; 73 : } 74 : 75 88 : double proj=0; mypack.clear(); 76 550 : for(unsigned j=0; j<vecs.size(); ++j) { 77 1848 : for(unsigned k=0; k<3; ++k) { 78 1386 : proj += vecs[j][k]*mypack.getAtomsDisplacementVector()[j][k]; 79 : } 80 462 : mypack.setAtomDerivatives( j, vecs[j] - comder ); 81 : } 82 88 : if( !mypack.updateComplete() ) mypack.updateDynamicLists(); 83 88 : return proj; 84 : } 85 : 86 : }