biod.pnpi.spb.ru
/
alexxy
/
gromacs.git
/ commitdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
| commitdiff |
tree
raw
|
patch
|
inline
| side by side (parent:
659677f
)
Pass ref_t to initialize_lambdas
author
Joe Jordan
<ejjordan12@gmail.com>
Sun, 11 Apr 2021 00:48:56 +0000
(
00:48
+0000)
committer
Andrey Alekseenko
<al42and@gmail.com>
Sun, 11 Apr 2021 00:48:56 +0000
(
00:48
+0000)
src/gromacs/mdrun/md.cpp
patch
|
blob
|
history
src/gromacs/mdrun/mimic.cpp
patch
|
blob
|
history
src/gromacs/mdrun/minimize.cpp
patch
|
blob
|
history
src/gromacs/mdrun/rerun.cpp
patch
|
blob
|
history
src/gromacs/mdtypes/state.cpp
patch
|
blob
|
history
src/gromacs/mdtypes/state.h
patch
|
blob
|
history
src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
patch
|
blob
|
history
src/gromacs/modularsimulator/freeenergyperturbationdata.h
patch
|
blob
|
history
src/gromacs/modularsimulator/simulatoralgorithm.cpp
patch
|
blob
|
history
diff --git
a/src/gromacs/mdrun/md.cpp
b/src/gromacs/mdrun/md.cpp
index 6e80825afef2ddc81dff8a3f7fda9f7d1f473587..a7487ff8cf79ef7f7dd245e223e74edfc8625f6e 100644
(file)
--- a/
src/gromacs/mdrun/md.cpp
+++ b/
src/gromacs/mdrun/md.cpp
@@
-257,7
+257,8
@@
void gmx::LegacySimulator::do_md()
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
- initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+ initialize_lambdas(
+ fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
Update upd(*ir, deform);
bool doSimulatedAnnealing = false;
{
Update upd(*ir, deform);
bool doSimulatedAnnealing = false;
{
diff --git
a/src/gromacs/mdrun/mimic.cpp
b/src/gromacs/mdrun/mimic.cpp
index ce23edc526ff30eb75d9c175c38209dd5ca70f1a..a434de0126b78830281683731e4fc81f90615ab6 100644
(file)
--- a/
src/gromacs/mdrun/mimic.cpp
+++ b/
src/gromacs/mdrun/mimic.cpp
@@
-228,7
+228,12
@@
void gmx::LegacySimulator::do_mimic()
nonConstGlobalTopology->intermolecularExclusionGroup = genQmmmIndices(top_global);
}
nonConstGlobalTopology->intermolecularExclusionGroup = genQmmmIndices(top_global);
}
- initialize_lambdas(fplog, *ir, MASTER(cr), &state_global->fep_state, state_global->lambda);
+ initialize_lambdas(fplog,
+ *ir,
+ gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc),
+ MASTER(cr),
+ &state_global->fep_state,
+ state_global->lambda);
const bool simulationsShareState = false;
gmx_mdoutf* outf = init_mdoutf(fplog,
const bool simulationsShareState = false;
gmx_mdoutf* outf = init_mdoutf(fplog,
diff --git
a/src/gromacs/mdrun/minimize.cpp
b/src/gromacs/mdrun/minimize.cpp
index 2083308bf01eaa385f6f9abe57e3e258706b90b3..c7da246c7cf15394c9c086e0386c08c74d119f51 100644
(file)
--- a/
src/gromacs/mdrun/minimize.cpp
+++ b/
src/gromacs/mdrun/minimize.cpp
@@
-394,7
+394,8
@@
static void init_em(FILE* fplog,
}
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
}
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
- initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+ initialize_lambdas(
+ fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
if (ir->eI == IntegrationAlgorithm::NM)
{
if (ir->eI == IntegrationAlgorithm::NM)
{
diff --git
a/src/gromacs/mdrun/rerun.cpp
b/src/gromacs/mdrun/rerun.cpp
index d3202d64fdd8f791a0527ad93bc33aeb601c14b4..ed4ede713415fbcef6297b86b9f560ee0e67f053 100644
(file)
--- a/
src/gromacs/mdrun/rerun.cpp
+++ b/
src/gromacs/mdrun/rerun.cpp
@@
-275,7
+275,8
@@
void gmx::LegacySimulator::do_rerun()
}
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
}
int* fep_state = MASTER(cr) ? &state_global->fep_state : nullptr;
gmx::ArrayRef<real> lambda = MASTER(cr) ? state_global->lambda : gmx::ArrayRef<real>();
- initialize_lambdas(fplog, *ir, MASTER(cr), fep_state, lambda);
+ initialize_lambdas(
+ fplog, *ir, gmx::arrayRefFromArray(ir->opts.ref_t, ir->opts.ngtc), MASTER(cr), fep_state, lambda);
const bool simulationsShareState = false;
gmx_mdoutf* outf = init_mdoutf(fplog,
nfile,
const bool simulationsShareState = false;
gmx_mdoutf* outf = init_mdoutf(fplog,
nfile,
diff --git
a/src/gromacs/mdtypes/state.cpp
b/src/gromacs/mdtypes/state.cpp
index 79a98f08455eb70017c5147bf25d4ca22ed555b6..4722245ff8570a13f4ec7b342b11258929953689 100644
(file)
--- a/
src/gromacs/mdtypes/state.cpp
+++ b/
src/gromacs/mdtypes/state.cpp
@@
-352,7
+352,12
@@
void printLambdaStateToLog(FILE* fplog, gmx::ArrayRef<const real> lambda, const
}
}
}
}
-void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* fep_state, gmx::ArrayRef<real> lambda)
+void initialize_lambdas(FILE* fplog,
+ const t_inputrec& ir,
+ gmx::ArrayRef<real> ref_t,
+ bool isMaster,
+ int* fep_state,
+ gmx::ArrayRef<real> lambda)
{
/* TODO: Clean up initialization of fep_state and lambda in
t_state. This function works, but could probably use a logic
{
/* TODO: Clean up initialization of fep_state and lambda in
t_state. This function works, but could probably use a logic
@@
-391,11
+396,11
@@
void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* f
if (ir.bSimTemp)
{
/* need to rescale control temperatures to match current state */
if (ir.bSimTemp)
{
/* need to rescale control temperatures to match current state */
- for (int i = 0; i <
ir.opts.ngtc
; i++)
+ for (int i = 0; i <
ref_t.ssize()
; i++)
{
{
- if (
ir.opts.
ref_t[i] > 0)
+ if (ref_t[i] > 0)
{
{
-
ir.opts.
ref_t[i] = ir.simtempvals->temperatures[fep->init_fep_state];
+ ref_t[i] = ir.simtempvals->temperatures[fep->init_fep_state];
}
}
}
}
}
}
diff --git
a/src/gromacs/mdtypes/state.h
b/src/gromacs/mdtypes/state.h
index 8998b75abcc86683c37af1d99ea4a0b774788405..faa51d21a45082bfeb4553cfaeaf46ea12010352 100644
(file)
--- a/
src/gromacs/mdtypes/state.h
+++ b/
src/gromacs/mdtypes/state.h
@@
-293,17
+293,6
@@
struct t_extmass
tensor Winvm; /* inverse pressure mass tensor, computed */
};
tensor Winvm; /* inverse pressure mass tensor, computed */
};
-
-typedef struct
-{
- real veta;
- double rscale;
- double vscale;
- double rvscale;
- double alpha;
- double* vscale_nhc;
-} t_vetavars;
-
#endif // DOXYGEN
//! Resizes the T- and P-coupling state variables
#endif // DOXYGEN
//! Resizes the T- and P-coupling state variables
@@
-376,6
+365,11
@@
void printLambdaStateToLog(FILE* fplog, gmx::ArrayRef<const real> lambda, bool i
* and lambda on master rank.
*
* Reports the initial lambda state to the log file. */
* and lambda on master rank.
*
* Reports the initial lambda state to the log file. */
-void initialize_lambdas(FILE* fplog, const t_inputrec& ir, bool isMaster, int* fep_state, gmx::ArrayRef<real> lambda);
+void initialize_lambdas(FILE* fplog,
+ const t_inputrec& ir,
+ gmx::ArrayRef<real> ref_t,
+ bool isMaster,
+ int* fep_state,
+ gmx::ArrayRef<real> lambda);
#endif
#endif
diff --git
a/src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
b/src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
index d780e69216eb38303ddadc503d75b77a04bfa11c..6671c20d85783871dfdec2700849df9f27f2a353 100644
(file)
--- a/
src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
+++ b/
src/gromacs/modularsimulator/freeenergyperturbationdata.cpp
@@
-60,8
+60,8
@@
namespace gmx
{
namespace gmx
{
-FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inputrec
*
inputrec, MDAtoms* mdAtoms) :
- element_(std::make_unique<Element>(this, inputrec
->
fepvals->delta_lambda)),
+FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inputrec
&
inputrec, MDAtoms* mdAtoms) :
+ element_(std::make_unique<Element>(this, inputrec
.
fepvals->delta_lambda)),
lambda_(),
currentFEPState_(0),
fplog_(fplog),
lambda_(),
currentFEPState_(0),
fplog_(fplog),
@@
-72,7
+72,12
@@
FreeEnergyPerturbationData::FreeEnergyPerturbationData(FILE* fplog, const t_inpu
// The legacy implementation only filled the lambda vector in state_global, which is only
// available on master. We have the lambda vector available everywhere, so we pass a `true`
// for isMaster on all ranks. See #3647.
// The legacy implementation only filled the lambda vector in state_global, which is only
// available on master. We have the lambda vector available everywhere, so we pass a `true`
// for isMaster on all ranks. See #3647.
- initialize_lambdas(fplog_, *inputrec_, true, ¤tFEPState_, lambda_);
+ initialize_lambdas(fplog_,
+ inputrec_,
+ gmx::arrayRefFromArray(inputrec_.opts.ref_t, inputrec_.opts.ngtc),
+ true,
+ ¤tFEPState_,
+ lambda_);
}
void FreeEnergyPerturbationData::Element::scheduleTask(Step step,
}
void FreeEnergyPerturbationData::Element::scheduleTask(Step step,
@@
-88,7
+93,7
@@
void FreeEnergyPerturbationData::Element::scheduleTask(Step step,
void FreeEnergyPerturbationData::updateLambdas(Step step)
{
// at beginning of step (if lambdas change...)
void FreeEnergyPerturbationData::updateLambdas(Step step)
{
// at beginning of step (if lambdas change...)
- lambda_ = currentLambdas(step, *(inputrec_
->
fepvals), currentFEPState_);
+ lambda_ = currentLambdas(step, *(inputrec_
.
fepvals), currentFEPState_);
updateMDAtoms();
}
updateMDAtoms();
}
@@
-202,7
+207,7
@@
void FreeEnergyPerturbationData::readCheckpointToTrxFrame(t_trxframe* trxFrame,
{
if (readCheckpointData)
{
{
if (readCheckpointData)
{
- FreeEnergyPerturbationData freeEnergyPerturbationData;
+ FreeEnergyPerturbationData freeEnergyPerturbationData
(nullptr, t_inputrec(), nullptr)
;
freeEnergyPerturbationData.doCheckpointData(&readCheckpointData.value());
trxFrame->lambda = freeEnergyPerturbationData.lambda_[FreeEnergyPerturbationCouplingType::Fep];
trxFrame->fep_state = freeEnergyPerturbationData.currentFEPState_;
freeEnergyPerturbationData.doCheckpointData(&readCheckpointData.value());
trxFrame->lambda = freeEnergyPerturbationData.lambda_[FreeEnergyPerturbationCouplingType::Fep];
trxFrame->fep_state = freeEnergyPerturbationData.currentFEPState_;
diff --git
a/src/gromacs/modularsimulator/freeenergyperturbationdata.h
b/src/gromacs/modularsimulator/freeenergyperturbationdata.h
index 7ef19d86f1e26fc21a8b4203dbb0081903356d97..98d9cfc80e5155a170783d0bf01288dfb9e583b7 100644
(file)
--- a/
src/gromacs/modularsimulator/freeenergyperturbationdata.h
+++ b/
src/gromacs/modularsimulator/freeenergyperturbationdata.h
@@
-77,7
+77,7
@@
class FreeEnergyPerturbationData final
{
public:
//! Constructor
{
public:
//! Constructor
- FreeEnergyPerturbationData(FILE* fplog, const t_inputrec
*
inputrec, MDAtoms* mdAtoms);
+ FreeEnergyPerturbationData(FILE* fplog, const t_inputrec
&
inputrec, MDAtoms* mdAtoms);
//! Get a view of the current lambda vector
ArrayRef<real> lambdaView();
//! Get a view of the current lambda vector
ArrayRef<real> lambdaView();
@@
-100,8
+100,6
@@
public:
static const std::string& checkpointID();
private:
static const std::string& checkpointID();
private:
- //! Default constructor - only used internally
- FreeEnergyPerturbationData() = default;
//! Update the lambda values
void updateLambdas(Step step);
//! Helper function to read from / write to CheckpointData
//! Update the lambda values
void updateLambdas(Step step);
//! Helper function to read from / write to CheckpointData
@@
-119,7
+117,7
@@
private:
//! Handles logging.
FILE* fplog_;
//! Contains user input mdp options.
//! Handles logging.
FILE* fplog_;
//! Contains user input mdp options.
- const t_inputrec
*
inputrec_;
+ const t_inputrec
&
inputrec_;
//! Atom parameters for this domain.
MDAtoms* mdAtoms_;
};
//! Atom parameters for this domain.
MDAtoms* mdAtoms_;
};
diff --git
a/src/gromacs/modularsimulator/simulatoralgorithm.cpp
b/src/gromacs/modularsimulator/simulatoralgorithm.cpp
index 423c3e463c2ce2f4a29c3254c58e3c0d48fdf71f..5c1106f1fdebc87732e7b052aa359b629db7d8af 100644
(file)
--- a/
src/gromacs/modularsimulator/simulatoralgorithm.cpp
+++ b/
src/gromacs/modularsimulator/simulatoralgorithm.cpp
@@
-407,7
+407,7
@@
ModularSimulatorAlgorithmBuilder::ModularSimulatorAlgorithmBuilder(
if (legacySimulatorData->inputrec->efep != FreeEnergyPerturbationType::No)
{
freeEnergyPerturbationData_ = std::make_unique<FreeEnergyPerturbationData>(
if (legacySimulatorData->inputrec->efep != FreeEnergyPerturbationType::No)
{
freeEnergyPerturbationData_ = std::make_unique<FreeEnergyPerturbationData>(
- legacySimulatorData->fplog, legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
+ legacySimulatorData->fplog,
*
legacySimulatorData->inputrec, legacySimulatorData->mdAtoms);
}
statePropagatorData_ = std::make_unique<StatePropagatorData>(
}
statePropagatorData_ = std::make_unique<StatePropagatorData>(