double MultiColvarBase::getCentralAtomDerivative( const unsigned& iatom, const unsigned& jcomp, const Vector& df ){ plumed_dbg_assert( atomsWithCatomDer.isActive(iatom) && jcomp<3 ); unsigned nder = 3*getNumberOfAtoms() + 9; return df[0]*getElementDerivative( (getCentralAtomElementIndex()+0)*nder + 3*iatom + jcomp ) + df[1]*getElementDerivative( (getCentralAtomElementIndex()+1)*nder + 3*iatom + jcomp ) + df[2]*getElementDerivative( (getCentralAtomElementIndex()+2)*nder + 3*iatom + jcomp ); }
void MultiColvarBase::mergeDerivatives( const unsigned& ider, const double& df ){ unsigned vstart=getNumberOfDerivatives()*ider; for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){ unsigned iatom=3*atoms_with_derivatives[i]; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); } unsigned nvir=3*getNumberOfAtoms(); for(unsigned j=0;j<9;++j){ accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++; } }
void MultiColvarBase::copyElementsToBridgedColvar( BridgedMultiColvarFunction* func ){ func->setElementValue( 0, getElementValue(0) ); for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){ unsigned n=atoms_with_derivatives[i], nx=3*n; func->atoms_with_derivatives.activate(n); func->addElementDerivative( nx+0, getElementDerivative(nx+0) ); func->addElementDerivative( nx+1, getElementDerivative(nx+1) ); func->addElementDerivative( nx+2, getElementDerivative(nx+2) ); } unsigned nvir=3*getNumberOfAtoms(); for(unsigned i=0;i<9;++i){ func->addElementDerivative( nvir, getElementDerivative(nvir) ); nvir++; } }
double OrientationSphere::compute(){ // Make sure derivatives for central atom are only calculated once VectorMultiColvar* vv = dynamic_cast<VectorMultiColvar*>( getBaseMultiColvar(0) ); vv->firstcall=true; weightHasDerivatives=true; // The weight has no derivatives really double sw, value=0, denom=0, dot, f_dot, dot_df, dfunc; Vector distance; getVectorForBaseTask(0, catom_orient ); for(unsigned i=1;i<getNAtoms();++i){ distance=getSeparation( getPositionOfCentralAtom(0), getPositionOfCentralAtom(i) ); sw = switchingFunction.calculateSqr( distance.modulo2(), dfunc ); if( sw>=getTolerance() ){ getVectorForBaseTask( i, this_orient ); // Calculate the dot product wrt to this position dot=0; for(unsigned k=0;k<catom_orient.size();++k) dot+=catom_orient[k]*this_orient[k]; f_dot = transformDotProduct( dot, dot_df ); // N.B. We are assuming here that the imaginary part of the dot product is zero for(unsigned k=0;k<catom_orient.size();++k){ this_orient[k]*=sw*dot_df; catom_der[k]=sw*dot_df*catom_orient[k]; } // Set the derivatives wrt of the numerator addOrientationDerivatives( 0, this_orient ); addOrientationDerivatives( i, catom_der ); addCentralAtomsDerivatives( 0, 0, f_dot*(-dfunc)*distance ); addCentralAtomsDerivatives( i, 0, f_dot*(dfunc)*distance ); addBoxDerivatives( f_dot*(-dfunc)*Tensor(distance,distance) ); value += sw*f_dot; // Set the derivatives wrt to the numerator addCentralAtomsDerivatives( 0, 1, (-dfunc)*distance ); addCentralAtomsDerivatives( i, 1, (dfunc)*distance ); addBoxDerivativesOfWeight( (-dfunc)*Tensor(distance,distance) ); denom += sw; } } // Now divide everything unsigned nder = getNumberOfDerivatives(); for(unsigned i=0;i<nder;++i){ setElementDerivative( i, getElementDerivative(i)/denom - (value*getElementDerivative(nder+i))/(denom*denom) ); setElementDerivative( nder + i, 0.0 ); } weightHasDerivatives=false; // Weight has no derivatives we just use the holder for weight to store some stuff return value / denom; }
void MultiColvarBase::quotientRule( const unsigned& uder, const unsigned& vder, const unsigned& iout ){ unsigned ustart=uder*getNumberOfDerivatives(); unsigned vstart=vder*getNumberOfDerivatives(); unsigned istart=iout*getNumberOfDerivatives(); double weight = getElementValue( vder ), pref = getElementValue( uder ) / (weight*weight); if( !doNotCalculateDerivatives() ){ for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){ unsigned n=3*atoms_with_derivatives[i], nx=n, ny=n+1, nz=n+2; setElementDerivative( istart + nx, getElementDerivative(ustart+nx) / weight - pref*getElementDerivative(vstart+nx) ); setElementDerivative( istart + ny, getElementDerivative(ustart+ny) / weight - pref*getElementDerivative(vstart+ny) ); setElementDerivative( istart + nz, getElementDerivative(ustart+nz) / weight - pref*getElementDerivative(vstart+nz) ); } unsigned vbase=3*getNumberOfAtoms(); for(unsigned i=0;i<9;++i){ setElementDerivative( istart + vbase + i, getElementDerivative(ustart+vbase+i) / weight - pref*getElementDerivative(vstart+vbase+i) ); } } thisval_wasset[iout]=false; setElementValue( iout, getElementValue(uder) / weight ); }
void BridgedMultiColvarFunction::mergeDerivatives( const unsigned& ider, const double& df ){ unsigned vstart=getNumberOfDerivatives()*ider; // Merge atom derivatives for(unsigned i=0;i<atoms_with_derivatives.getNumberActive();++i){ unsigned iatom=3*atoms_with_derivatives[i]; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); iatom++; accumulateDerivative( iatom, df*getElementDerivative(vstart+iatom) ); } // Merge virial derivatives unsigned nvir=3*mycolv->getNumberOfAtoms(); for(unsigned j=0;j<9;++j){ accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++; } // Merge local atom derivatives for(unsigned j=0;j<getNumberOfAtoms();++j){ accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++; accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++; accumulateDerivative( nvir, df*getElementDerivative(vstart+nvir) ); nvir++; } plumed_dbg_assert( nvir==getNumberOfDerivatives() ); }