45#include "flucq/FluctuatingChargeObjectiveFunction.hpp"
53 FluctuatingChargeObjectiveFunction::FluctuatingChargeObjectiveFunction(
57 forceMan_(forceMan), fqConstraints_(fqConstraints), thermo(info) {}
60 const DynamicVector<RealType>& x) {
62 forceMan_->calcForces();
63 fqConstraints_->applyConstraints();
64 return thermo.getPotential();
68 DynamicVector<RealType>& grad,
const DynamicVector<RealType>& x) {
70 forceMan_->calcForces();
71 fqConstraints_->applyConstraints();
76 DynamicVector<RealType>& grad,
const DynamicVector<RealType>& x) {
78 forceMan_->calcForces();
79 fqConstraints_->applyConstraints();
81 return thermo.getPotential();
84 void FluctuatingChargeObjectiveFunction::setCoor(
85 const DynamicVector<RealType>& x)
const {
86 SimInfo::MoleculeIterator i;
87 Molecule::FluctuatingChargeIterator j;
95 index = displacements_[myrank_];
102 for (atom = mol->beginFluctuatingCharge(j); atom != NULL;
103 atom = mol->nextFluctuatingCharge(j)) {
104 atom->setFlucQPos(x[index++]);
109 void FluctuatingChargeObjectiveFunction::getGrad(
110 DynamicVector<RealType>& grad) {
111 SimInfo::MoleculeIterator i;
112 Molecule::FluctuatingChargeIterator j;
119 index = displacements_[myrank_];
126 for (atom = mol->beginFluctuatingCharge(j); atom != NULL;
127 atom = mol->nextFluctuatingCharge(j)) {
128 grad[index++] = -atom->getFlucQFrc();
133 MPI_Allreduce(MPI_IN_PLACE, &grad[0], nFlucQ_, MPI_REALTYPE, MPI_SUM,
138 DynamicVector<RealType>
139 FluctuatingChargeObjectiveFunction::setInitialCoords() {
141 MPI_Comm_size(MPI_COMM_WORLD, &nproc_);
142 MPI_Comm_rank(MPI_COMM_WORLD, &myrank_);
143 std::vector<int> flucqOnProc_(nproc_, 0);
145 displacements_.clear();
146 displacements_.resize(nproc_, 0);
149 SimInfo::MoleculeIterator i;
150 Molecule::FluctuatingChargeIterator j;
158 for (atom = mol->beginFluctuatingCharge(j); atom != NULL;
159 atom = mol->nextFluctuatingCharge(j)) {
165 MPI_Allgather(&nFlucQ_, 1, MPI_INT, &flucqOnProc_[0], 1, MPI_INT,
169 for (
int iproc = 0; iproc < nproc_; iproc++) {
170 nFlucQ_ += flucqOnProc_[iproc];
173 displacements_[0] = 0;
174 for (
int iproc = 1; iproc < nproc_; iproc++) {
175 displacements_[iproc] =
176 displacements_[iproc - 1] + flucqOnProc_[iproc - 1];
180 DynamicVector<RealType> initCoords(nFlucQ_);
184 index = displacements_[myrank_];
191 for (atom = mol->beginFluctuatingCharge(j); atom != NULL;
192 atom = mol->nextFluctuatingCharge(j)) {
193 initCoords[index++] = atom->getFlucQPos();
198 MPI_Allgatherv(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &initCoords[0],
199 &flucqOnProc_[0], &displacements_[0], MPI_REALTYPE,
virtual RealType value(const DynamicVector< RealType > &x)
method to overload to compute the objective function value in x
virtual void gradient(DynamicVector< RealType > &grad, const DynamicVector< RealType > &x)
method to overload to compute grad_f, the first derivative of
virtual RealType valueAndGradient(DynamicVector< RealType > &grad, const DynamicVector< RealType > &x)
method to overload to compute grad_f, the first derivative
ForceManager is responsible for calculating both the short range (bonded) interactions and long range...
One of the heavy-weight classes of OpenMD, SimInfo maintains objects and variables relating to the cu...
Molecule * beginMolecule(MoleculeIterator &i)
Returns the first molecule in this SimInfo and intialize the iterator.
Molecule * nextMolecule(MoleculeIterator &i)
Returns the next avaliable Molecule based on the iterator.
SnapshotManager * getSnapshotManager()
Returns the snapshot manager.
This basic Periodic Table class was originally taken from the data.cpp file in OpenBabel.