diff --git a/src/multicolvar/InPlaneDistances.cpp b/src/multicolvar/InPlaneDistances.cpp new file mode 100644 index 0000000000000000000000000000000000000000..bf5ef362596011f1d728cd7e79f8d6fcd43ce456 --- /dev/null +++ b/src/multicolvar/InPlaneDistances.cpp @@ -0,0 +1,133 @@ +/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + Copyright (c) 2013,2014 The plumed team + (see the PEOPLE file at the root of the distribution for a list of names) + + See http://www.plumed-code.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 "MultiColvar.h" +#include "tools/SwitchingFunction.h" +#include "core/ActionRegister.h" +#include "vesselbase/LessThan.h" +#include "vesselbase/Between.h" +#include "tools/Angle.h" + +#include <string> +#include <cmath> + +using namespace std; + +namespace PLMD{ +namespace multicolvar{ + +//+PLUMEDOC MCOLVAR INPLANEDISTANCES +/* +Calculate distances in the plane perpendicular to an axis + +\par Examples + +*/ +//+ENDPLUMEDOC + +class InPlaneDistances : public MultiColvar { +public: + static void registerKeywords( Keywords& keys ); + InPlaneDistances(const ActionOptions&); +// active methods: + virtual double compute(); + bool isPeriodic(){ return false; } + Vector getCentralAtom(); +}; + +PLUMED_REGISTER_ACTION(InPlaneDistances,"INPLANEDISTANCES") + +void InPlaneDistances::registerKeywords( Keywords& keys ){ + MultiColvar::registerKeywords( keys ); + keys.use("MEAN"); keys.use("MIN"); keys.use("MAX"); keys.use("LESS_THAN"); + keys.use("MORE_THAN"); keys.use("BETWEEN"); keys.use("HISTOGRAM"); keys.use("MOMENTS"); + keys.add("atoms","VECTORSTART","The first atom position that is used to define the normal to the plane of interest"); + keys.add("atoms","VECTOREND","The second atom position that is used to defin the normal to the plane of interest"); + keys.add("atoms-2","GROUP","The set of atoms for which you wish to calculate the in plane distance "); +} + +InPlaneDistances::InPlaneDistances(const ActionOptions&ao): +PLUMED_MULTICOLVAR_INIT(ao) +{ + // Read in the atoms + std::vector<AtomNumber> all_atoms; + readThreeGroups("GROUP","VECTORSTART","VECTOREND",false, all_atoms); + if( all_atoms.size()>0 ) ActionAtomistic::requestAtoms( all_atoms ); + + // Check atoms are OK + if( getFullNumberOfTasks()!=getNumberOfAtoms()-2 ) error("you should specify one atom for VECTORSTART and one atom for VECTOREND only"); + + // Setup the multicolvar base + setupMultiColvarBase(); readVesselKeywords(); + // And check everything has been read in correctly + checkRead(); + + // Now check if we can use link cells + bool use_link=false; double rcut; + if( getNumberOfVessels()>0 ){ + vesselbase::LessThan* lt=dynamic_cast<vesselbase::LessThan*>( getPntrToVessel(0) ); + if( lt ){ + use_link=true; rcut=lt->getCutoff(); + } else { + vesselbase::Between* bt=dynamic_cast<vesselbase::Between*>( getPntrToVessel(0) ); + if( bt ) use_link=true; rcut=bt->getCutoff(); + } + if( use_link ){ + for(unsigned i=1;i<getNumberOfVessels();++i){ + vesselbase::LessThan* lt2=dynamic_cast<vesselbase::LessThan*>( getPntrToVessel(i) ); + vesselbase::Between* bt=dynamic_cast<vesselbase::Between*>( getPntrToVessel(i) ); + if( lt2 ){ + double tcut=lt2->getCutoff(); + if( tcut>rcut ) rcut=tcut; + } else if( bt ){ + double tcut=bt->getCutoff(); + if( tcut>rcut ) rcut=tcut; + } else { + use_link=false; + } + } + } + if( use_link ) setLinkCellCutoff( rcut ); + } +} + +double InPlaneDistances::compute(){ + Vector normal=getSeparation( getPosition(1), getPosition(2) ); + Vector dir=getSeparation( getPosition(1), getPosition(0) ); + PLMD::Angle a; Vector ddij, ddik; double angle=a.compute(normal,dir,ddij,ddik); + double sangle=sin(angle), cangle=cos(angle); + double dd=dir.modulo(), invdd=1.0/dd, val=dd*sangle; + + addAtomsDerivatives( 0, dd*cangle*ddik + sangle*invdd*dir ); + addAtomsDerivatives( 1, -dd*cangle*(ddik+ddij) - sangle*invdd*dir ); + addAtomsDerivatives( 2, dd*cangle*ddij ); + addBoxDerivatives( -dd*cangle*(Tensor(normal,ddij)+Tensor(dir,ddik)) - sangle*invdd*Tensor(dir,dir) ); + + return val; +} + +Vector InPlaneDistances::getCentralAtom(){ + addCentralAtomDerivatives( 0, Tensor::identity() ); + return getPosition(0); +} + +} +} diff --git a/src/reference/DotProductDistance.cpp b/src/reference/DotProductDistance.cpp new file mode 100644 index 0000000000000000000000000000000000000000..7dd4d8f8697eed4063ecb9b3396c85862433fb3c --- /dev/null +++ b/src/reference/DotProductDistance.cpp @@ -0,0 +1,54 @@ +/* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + Copyright (c) 2012-2014 The plumed team + (see the PEOPLE file at the root of the distribution for a list of names) + + See http://www.plumed-code.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 "MetricRegister.h" +#include "ArgumentOnlyDistance.h" +#include "core/Value.h" + +namespace PLMD { + +class DotProductDistance : public ArgumentOnlyDistance { +public: + DotProductDistance( const ReferenceConfigurationOptions& ro ); + void read( const PDB& ); + double calc( const std::vector<Value*>& vals, const std::vector<double>& arg, const bool& squared ); +}; + +PLUMED_REGISTER_METRIC(DotProductDistance,"DOTPRODUCT") + +DotProductDistance::DotProductDistance( const ReferenceConfigurationOptions& ro ): +ReferenceConfiguration(ro), +ArgumentOnlyDistance(ro) +{ +} + +void DotProductDistance::read( const PDB& pdb ){ + readArgumentsFromPDB( pdb ); +} + +double DotProductDistance::calc( const std::vector<Value*>& vals, const std::vector<double>& arg, const bool& squared ){ + double dot=0.0; + for (unsigned long i=0; i<vals.size(); ++i) dot+=getReferenceArgument(i)*arg[i]; + return -log(dot); +} + + +}