The previous implementation would lead to slight inaccuracies when
virtual positions were calculated before shifting positions in
pressure coupling algorithms. This change moves the calculation to the
top of the loop, right before DD and force calculation.
This also removes an unnecessary call to the virtual site construction
in rerun. That call was placed after force calculation and trajectory
writing, so had no influence since coordinates for the next step are
read from file. Since virtual sites are constructed when reading
coordinates from file, the results should have been correct anyway, the
call just had no influence.
This also moves the call to the virtual site construction in MiMiC.
Similar to the rerun case, the current call happened after virtual
sites are used. This call has been moved to right after the coordinates
are received. This is likely to have lead to wrong results when using
MiMiC with virtual sites.
Refs #3866
+Fixed slight inaccuracies when using virtual sites with pressure coupling
+"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
+
+Virtual sites were reconstructed after the system was propagated, but before
+scaling due to pressure coupling. For virtual site types which are not a linear
+combination of other atoms, this is not completely correct. Since the scaling
+due to pressure coupling is very small in healthy simulations, the resulting
+inaccuracies are expected to have been extremely minor, and in most cases
+undetectable.
+
+:issue:`3866`
+
.. Note to developers!
Please use """"""" to underline the individual entries for fixed issues in the subfolders,
otherwise the formatting on the webpage is messed up.
.. Note to developers!
Please use """"""" to underline the individual entries for fixed issues in the subfolders,
otherwise the formatting on the webpage is messed up.
state->box,
state->lambda[efptBONDED]);
}
state->box,
state->lambda[efptBONDED]);
}
- if (vsite)
- {
- /* Construct the virtual sites for the initial configuration */
- vsite->construct(state->x, ir->delta_t, {}, state->box);
- }
}
if (ir->efep != efepNO)
}
if (ir->efep != efepNO)
stateGpu->waitCoordinatesReadyOnHost(AtomLocality::Local);
}
stateGpu->waitCoordinatesReadyOnHost(AtomLocality::Local);
}
+ if (vsite != nullptr)
+ {
+ // Virtual sites need to be updated before domain decomposition and forces are calculated
+ wallcycle_start(wcycle, ewcVSITECONSTR);
+ vsite->construct(state->x, ir->delta_t, state->v, state->box);
+ wallcycle_stop(wcycle, ewcVSITECONSTR);
+ }
+
if (bNS && !(bFirstStep && ir->bContinuation))
{
bMasterState = FALSE;
if (bNS && !(bFirstStep && ir->bContinuation))
{
bMasterState = FALSE;
enerd->term[F_DVDL_CONSTR] += dvdl_constr;
}
enerd->term[F_DVDL_CONSTR] += dvdl_constr;
}
- if (vsite != nullptr)
- {
- wallcycle_start(wcycle, ewcVSITECONSTR);
- vsite->construct(state->x, ir->delta_t, state->v, state->box);
- wallcycle_stop(wcycle, ewcVSITECONSTR);
- }
-
/* ############## IF NOT VV, Calculate globals HERE ############ */
/* With Leap-Frog we can skip compute_globals at
* non-communication steps, but we need to calculate
/* ############## IF NOT VV, Calculate globals HERE ############ */
/* With Leap-Frog we can skip compute_globals at
* non-communication steps, but we need to calculate
"decomposition, "
"use a single rank");
}
"decomposition, "
"use a single rank");
}
+ if (constructVsites)
+ {
+ wallcycle_start(wcycle, ewcVSITECONSTR);
+ vsite->construct(state->x, ir->delta_t, state->v, state->box);
+ wallcycle_stop(wcycle, ewcVSITECONSTR);
+ }
stopHandler->setSignal();
stopHandler->setSignal();
- if (vsite != nullptr)
- {
- wallcycle_start(wcycle, ewcVSITECONSTR);
- vsite->construct(state->x, ir->delta_t, state->v, state->box);
- wallcycle_stop(wcycle, ewcVSITECONSTR);
- }
-
{
const bool doInterSimSignal = false;
const bool doIntraSimSignal = true;
{
const bool doInterSimSignal = false;
const bool doIntraSimSignal = true;
stopHandler->setSignal();
stopHandler->setSignal();
- if (vsite != nullptr)
- {
- wallcycle_start(wcycle, ewcVSITECONSTR);
- vsite->construct(state->x, ir->delta_t, state->v, state->box);
- wallcycle_stop(wcycle, ewcVSITECONSTR);
- }
-
{
const bool doInterSimSignal = false;
const bool doIntraSimSignal = true;
{
const bool doInterSimSignal = false;
const bool doIntraSimSignal = true;