Skip to content

Commit

Permalink
Merge pull request #281 from OpenBioSim/feature_membrane_barostat
Browse files Browse the repository at this point in the history 
Add support for OpenMM MonteCarloMembraneBarostat
  • Loading branch information
@lohedges
lohedges authored Feb 4, 2025
2 parents b3a566b + b29e168 commit a5278e9
Show file tree
Hide file tree
Showing 14 changed files with 284 additions and 44 deletions.
45 changes: 37 additions & 8 deletions corelib/src/libs/SireMove/openmmfrenergydt.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -129,7 +129,7 @@ QDataStream &operator<<(QDataStream &ds, const OpenMMFrEnergyDT &velver)

sds << velver.frequent_save_velocities << velver.molgroup << velver.solutegroup << velver.CutoffType
<< velver.cutoff_distance << velver.field_dielectric << velver.Andersen_flag << velver.Andersen_frequency
<< velver.MCBarostat_flag << velver.MCBarostat_frequency << velver.ConstraintType << velver.Pressure
<< velver.MCBarostat_flag << velver.MCBarostat_membrane_flag << velver.MCBarostat_frequency << velver.ConstraintType << velver.Pressure
<< velver.Temperature << velver.platform_type << velver.Restraint_flag << velver.CMMremoval_frequency
<< velver.energy_frequency << velver.device_index << velver.Alchemical_value << velver.coulomb_power
<< velver.shift_delta << velver.delta_alchemical << velver.buffer_coords << velver.gradients
Expand All @@ -151,7 +151,7 @@ QDataStream &operator>>(QDataStream &ds, OpenMMFrEnergyDT &velver)

sds >> velver.frequent_save_velocities >> velver.molgroup >> velver.solutegroup >> velver.CutoffType >>
velver.cutoff_distance >> velver.field_dielectric >> velver.Andersen_flag >> velver.Andersen_frequency >>
velver.MCBarostat_flag >> velver.MCBarostat_frequency >> velver.ConstraintType >> velver.Pressure >>
velver.MCBarostat_flag >> velver.MCBarostat_membrane_flag >> velver.MCBarostat_frequency >> velver.ConstraintType >> velver.Pressure >>
velver.Temperature >> velver.platform_type >> velver.Restraint_flag >> velver.CMMremoval_frequency >>
velver.energy_frequency >> velver.device_index >> velver.Alchemical_value >> velver.coulomb_power >>
velver.shift_delta >> velver.delta_alchemical >> velver.buffer_coords >> velver.gradients >>
Expand All @@ -175,7 +175,7 @@ OpenMMFrEnergyDT::OpenMMFrEnergyDT(bool frequent_save)
: ConcreteProperty<OpenMMFrEnergyDT, Integrator>(), frequent_save_velocities(frequent_save),
molgroup(MoleculeGroup()), solutegroup(MoleculeGroup()), openmm_system(0), isInitialised(false),
CutoffType("nocutoff"), cutoff_distance(1.0 * nanometer), field_dielectric(78.3), Andersen_flag(false),
Andersen_frequency(90.0), MCBarostat_flag(false), MCBarostat_frequency(25), ConstraintType("none"),
Andersen_frequency(90.0), MCBarostat_flag(false), MCBarostat_membrane_flag(false), MCBarostat_frequency(25), ConstraintType("none"),
Pressure(1.0 * bar), Temperature(300.0 * kelvin), platform_type("Reference"), Restraint_flag(false),
CMMremoval_frequency(0), energy_frequency(100), device_index("0"), Alchemical_value(0.5), coulomb_power(0),
shift_delta(2.0), delta_alchemical(0.001), buffer_coords(false), gradients()
Expand All @@ -188,7 +188,7 @@ OpenMMFrEnergyDT::OpenMMFrEnergyDT(const MoleculeGroup &molecule_group, const Mo
: ConcreteProperty<OpenMMFrEnergyDT, Integrator>(), frequent_save_velocities(frequent_save),
molgroup(molecule_group), solutegroup(solute_group), openmm_system(0), isInitialised(false),
CutoffType("nocutoff"), cutoff_distance(1.0 * nanometer), field_dielectric(78.3), Andersen_flag(false),
Andersen_frequency(90.0), MCBarostat_flag(false), MCBarostat_frequency(25), ConstraintType("none"),
Andersen_frequency(90.0), MCBarostat_flag(false), MCBarostat_membrane_flag(false), MCBarostat_frequency(25), ConstraintType("none"),
Pressure(1.0 * bar), Temperature(300.0 * kelvin), platform_type("Reference"), Restraint_flag(false),
CMMremoval_frequency(0), energy_frequency(100), device_index("0"), Alchemical_value(0.5), coulomb_power(0),
shift_delta(2.0), delta_alchemical(0.001), buffer_coords(false), gradients()
Expand All @@ -201,7 +201,7 @@ OpenMMFrEnergyDT::OpenMMFrEnergyDT(const OpenMMFrEnergyDT &other)
molgroup(other.molgroup), solutegroup(other.solutegroup), openmm_system(other.openmm_system),
isInitialised(other.isInitialised), CutoffType(other.CutoffType), cutoff_distance(other.cutoff_distance),
field_dielectric(other.field_dielectric), Andersen_flag(other.Andersen_flag),
Andersen_frequency(other.Andersen_frequency), MCBarostat_flag(other.MCBarostat_flag),
Andersen_frequency(other.Andersen_frequency), MCBarostat_flag(other.MCBarostat_flag), MCBarostat_membrane_flag(other.MCBarostat_membrane_flag),
MCBarostat_frequency(other.MCBarostat_frequency), ConstraintType(other.ConstraintType), Pressure(other.Pressure),
Temperature(other.Temperature), platform_type(other.platform_type), Restraint_flag(other.Restraint_flag),
CMMremoval_frequency(other.CMMremoval_frequency), energy_frequency(other.energy_frequency),
Expand Down Expand Up @@ -232,6 +232,7 @@ OpenMMFrEnergyDT &OpenMMFrEnergyDT::operator=(const OpenMMFrEnergyDT &other)
Andersen_flag = other.Andersen_flag;
Andersen_frequency = other.Andersen_frequency;
MCBarostat_flag = other.MCBarostat_flag;
MCBarostat_membrane_flag = other.MCBarostat_membrane_flag;
MCBarostat_frequency = other.MCBarostat_frequency;
ConstraintType = other.ConstraintType;
Pressure = other.Pressure;
Expand Down Expand Up @@ -521,16 +522,32 @@ void OpenMMFrEnergyDT::initialise()
const double converted_Temperature = convertTo(Temperature.value(), kelvin);
const double converted_Pressure = convertTo(Pressure.value(), bar);

OpenMM::MonteCarloBarostat *barostat =
new OpenMM::MonteCarloBarostat(converted_Pressure, converted_Temperature, MCBarostat_frequency);
system_openmm->addForce(barostat);
if (MCBarostat_membrane_flag)
{
// Simple options for now: zero surface tension, XY isotropic, Z free
const double surface_Tension = 0;
OpenMM::MonteCarloMembraneBarostat::XYMode xymode = OpenMM::MonteCarloMembraneBarostat::XYIsotropic;
OpenMM::MonteCarloMembraneBarostat::ZMode zmode = OpenMM::MonteCarloMembraneBarostat::ZFree;
OpenMM::MonteCarloMembraneBarostat * barostat = new OpenMM::MonteCarloMembraneBarostat(converted_Pressure, surface_Tension, converted_Temperature, xymode, zmode, MCBarostat_frequency);
system_openmm->addForce(barostat);
}
else
{
OpenMM::MonteCarloBarostat *barostat =
new OpenMM::MonteCarloBarostat(converted_Pressure, converted_Temperature, MCBarostat_frequency);
system_openmm->addForce(barostat);
}

if (Debug)
{
qDebug() << "\nMonte Carlo Barostat set\n";
qDebug() << "Temperature = " << converted_Temperature << " K\n";
qDebug() << "Pressure = " << converted_Pressure << " bar\n";
qDebug() << "Frequency every " << MCBarostat_frequency << " steps\n";
if (MCBarostat_membrane_flag)
{
qDebug() << "Membrane barostat, surface tension 0, XY isotropic, Z free\n";
}
}
}

Expand Down Expand Up @@ -1795,6 +1812,18 @@ bool OpenMMFrEnergyDT::getMCBarostat(void)
return MCBarostat_flag;
}

/** Set Monte Carlo membrane Barostat on/off */
void OpenMMFrEnergyDT::setMCBarostat_membrane(bool MCBarostat_membrane)
{
MCBarostat_membrane_flag = MCBarostat_membrane;
}

bool OpenMMFrEnergyDT::getMCBarostat_membrane(void)
{

return MCBarostat_membrane_flag;
}

/** Get the Monte Carlo Barostat frequency in time speps */
int OpenMMFrEnergyDT::getMCBarostat_frequency(void)
{
Expand Down
4 changes: 4 additions & 0 deletions corelib/src/libs/SireMove/openmmfrenergydt.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -114,6 +114,9 @@ namespace SireMove
void setMCBarostat_frequency(int);
int getMCBarostat_frequency(void);

bool getMCBarostat_membrane(void);
void setMCBarostat_membrane(bool);

QString getConstraintType(void);
void setConstraintType(QString);

Expand Down Expand Up @@ -180,6 +183,7 @@ namespace SireMove
double Andersen_frequency;

bool MCBarostat_flag;
bool MCBarostat_membrane_flag;
int MCBarostat_frequency;

QString ConstraintType;
Expand Down
53 changes: 42 additions & 11 deletions corelib/src/libs/SireMove/openmmfrenergyst.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -130,7 +130,7 @@ QDataStream &operator<<(QDataStream &ds, const OpenMMFrEnergyST &velver)
sds << velver.frequent_save_velocities << velver.molgroup << velver.solute << velver.solutehard
<< velver.solutetodummy << velver.solutefromdummy << velver.combiningRules << velver.CutoffType
<< velver.cutoff_distance << velver.field_dielectric << velver.Andersen_flag << velver.Andersen_frequency
<< velver.MCBarostat_flag << velver.MCBarostat_frequency << velver.ConstraintType << velver.Pressure
<< velver.MCBarostat_flag << velver.MCBarostat_membrane_flag << velver.MCBarostat_frequency << velver.ConstraintType << velver.Pressure
<< velver.Temperature << velver.platform_type << velver.Restraint_flag << velver.CMMremoval_frequency
<< velver.buffer_frequency << velver.energy_frequency << velver.device_index << velver.precision
<< velver.Alchemical_value << velver.coulomb_power << velver.shift_delta << velver.delta_alchemical
Expand All @@ -152,7 +152,7 @@ QDataStream &operator>>(QDataStream &ds, OpenMMFrEnergyST &velver)
sds >> velver.frequent_save_velocities >> velver.molgroup >> velver.solute >> velver.solutehard >>
velver.solutetodummy >> velver.solutefromdummy >> velver.combiningRules >> velver.CutoffType >>
velver.cutoff_distance >> velver.field_dielectric >> velver.Andersen_flag >> velver.Andersen_frequency >>
velver.MCBarostat_flag >> velver.MCBarostat_frequency >> velver.ConstraintType >> velver.Pressure >>
velver.MCBarostat_flag >> velver.MCBarostat_membrane_flag >> velver.MCBarostat_frequency >> velver.ConstraintType >> velver.Pressure >>
velver.Temperature >> velver.platform_type >> velver.Restraint_flag >> velver.CMMremoval_frequency >>
velver.buffer_frequency >> velver.energy_frequency >> velver.device_index >> velver.precision >>
velver.Alchemical_value >> velver.coulomb_power >> velver.shift_delta >> velver.delta_alchemical >>
Expand Down Expand Up @@ -180,7 +180,7 @@ OpenMMFrEnergyST::OpenMMFrEnergyST(bool frequent_save)
solutefromdummy(MoleculeGroup()), openmm_system(0), openmm_context(0), isSystemInitialised(false),
isContextInitialised(false), combiningRules("arithmetic"), CutoffType("nocutoff"),
cutoff_distance(1.0 * nanometer), field_dielectric(78.3), Andersen_flag(false), Andersen_frequency(90.0),
MCBarostat_flag(false), MCBarostat_frequency(25), ConstraintType("none"), Pressure(1.0 * bar),
MCBarostat_flag(false), MCBarostat_membrane_flag(false), MCBarostat_frequency(25), ConstraintType("none"), Pressure(1.0 * bar),
Temperature(300.0 * kelvin), platform_type("Reference"), Restraint_flag(false), CMMremoval_frequency(0),
buffer_frequency(0), energy_frequency(100), device_index("0"), precision("single"), Alchemical_value(0.5),
coulomb_power(0), shift_delta(2.0), delta_alchemical(0.001), alchemical_array(), finite_diff_gradients(),
Expand All @@ -199,7 +199,7 @@ OpenMMFrEnergyST::OpenMMFrEnergyST(const MoleculeGroup &molecule_group, const Mo
solutefromdummy(solute_fromdummy), openmm_system(0), openmm_context(0), isSystemInitialised(false),
isContextInitialised(false), combiningRules("arithmetic"), CutoffType("nocutoff"),
cutoff_distance(1.0 * nanometer), field_dielectric(78.3), Andersen_flag(false), Andersen_frequency(90.0),
MCBarostat_flag(false), MCBarostat_frequency(25), ConstraintType("none"), Pressure(1.0 * bar),
MCBarostat_flag(false), MCBarostat_membrane_flag(false), MCBarostat_frequency(25), ConstraintType("none"), Pressure(1.0 * bar),
Temperature(300.0 * kelvin), platform_type("Reference"), Restraint_flag(false), CMMremoval_frequency(0),
buffer_frequency(0), energy_frequency(100), device_index("0"), precision("single"), Alchemical_value(0.5),
coulomb_power(0), shift_delta(2.0), delta_alchemical(0.001), alchemical_array(), finite_diff_gradients(),
Expand All @@ -217,7 +217,7 @@ OpenMMFrEnergyST::OpenMMFrEnergyST(const OpenMMFrEnergyST &other)
isSystemInitialised(other.isSystemInitialised), isContextInitialised(other.isContextInitialised),
combiningRules(other.combiningRules), CutoffType(other.CutoffType), cutoff_distance(other.cutoff_distance),
field_dielectric(other.field_dielectric), Andersen_flag(other.Andersen_flag),
Andersen_frequency(other.Andersen_frequency), MCBarostat_flag(other.MCBarostat_flag),
Andersen_frequency(other.Andersen_frequency), MCBarostat_flag(other.MCBarostat_flag), MCBarostat_membrane_flag(other.MCBarostat_membrane_flag),
MCBarostat_frequency(other.MCBarostat_frequency), ConstraintType(other.ConstraintType), Pressure(other.Pressure),
Temperature(other.Temperature), platform_type(other.platform_type), Restraint_flag(other.Restraint_flag),
CMMremoval_frequency(other.CMMremoval_frequency), buffer_frequency(other.buffer_frequency),
Expand Down Expand Up @@ -259,6 +259,7 @@ OpenMMFrEnergyST &OpenMMFrEnergyST::operator=(const OpenMMFrEnergyST &other)
Andersen_flag = other.Andersen_flag;
Andersen_frequency = other.Andersen_frequency;
MCBarostat_flag = other.MCBarostat_flag;
MCBarostat_membrane_flag = other.MCBarostat_membrane_flag;
MCBarostat_frequency = other.MCBarostat_frequency;
ConstraintType = other.ConstraintType;
Pressure = other.Pressure;
Expand Down Expand Up @@ -302,7 +303,7 @@ bool OpenMMFrEnergyST::operator==(const OpenMMFrEnergyST &other) const
combiningRules == other.combiningRules and CutoffType == other.CutoffType and
cutoff_distance == other.cutoff_distance and field_dielectric == other.field_dielectric and
Andersen_flag == other.Andersen_flag and Andersen_frequency == other.Andersen_frequency and
MCBarostat_flag == other.MCBarostat_flag and MCBarostat_frequency == other.MCBarostat_frequency and
MCBarostat_flag == other.MCBarostat_flag and MCBarostat_membrane_flag == other.MCBarostat_membrane_flag and MCBarostat_frequency == other.MCBarostat_frequency and
ConstraintType == other.ConstraintType and Pressure == other.Pressure and
Temperature == other.Temperature and platform_type == other.platform_type and
Restraint_flag == other.Restraint_flag and CMMremoval_frequency == other.CMMremoval_frequency and
Expand Down Expand Up @@ -1241,20 +1242,40 @@ void OpenMMFrEnergyST::initialise()
const double converted_Temperature = convertTo(Temperature.value(), kelvin);
const double converted_Pressure = convertTo(Pressure.value(), bar);

OpenMM::MonteCarloBarostat *barostat =
new OpenMM::MonteCarloBarostat(converted_Pressure, converted_Temperature, MCBarostat_frequency);
if (MCBarostat_membrane_flag == true)
{
// Simple options for now: zero surface tension, XY isotropic, Z free
const double surface_Tension = 0;
OpenMM::MonteCarloMembraneBarostat::XYMode xymode = OpenMM::MonteCarloMembraneBarostat::XYIsotropic;
OpenMM::MonteCarloMembraneBarostat::ZMode zmode = OpenMM::MonteCarloMembraneBarostat::ZFree;
OpenMM::MonteCarloMembraneBarostat * barostat = new OpenMM::MonteCarloMembraneBarostat(converted_Pressure, surface_Tension, converted_Temperature, xymode, zmode, MCBarostat_frequency);

// Set The random seed
barostat->setRandomNumberSeed(random_seed);
//Set The random seed
barostat->setRandomNumberSeed(random_seed);

system_openmm->addForce(barostat);
}
else
{
OpenMM::MonteCarloBarostat *barostat =
new OpenMM::MonteCarloBarostat(converted_Pressure, converted_Temperature, MCBarostat_frequency);

system_openmm->addForce(barostat);
// Set The random seed
barostat->setRandomNumberSeed(random_seed);

system_openmm->addForce(barostat);
}

if (Debug)
{
qDebug() << "\nMonte Carlo Barostat set\n";
qDebug() << "Temperature = " << converted_Temperature << " K\n";
qDebug() << "Pressure = " << converted_Pressure << " bar\n";
qDebug() << "Frequency every " << MCBarostat_frequency << " steps\n";
if (MCBarostat_membrane_flag)
{
qDebug() << "Membrane barostat, surface tension 0, XY isotropic, Z free\n";
}
}
}
/*******************************************************BONDED
Expand Down Expand Up @@ -4351,6 +4372,16 @@ bool OpenMMFrEnergyST::getMCBarostat(void)
return MCBarostat_flag;
}

void OpenMMFrEnergyST::setMCBarostatMembrane(bool MCBarostat_membrane)
{
MCBarostat_membrane_flag = MCBarostat_membrane;
}

bool OpenMMFrEnergyST::getMCBarostatMembrane(void)
{
return MCBarostat_membrane_flag;
}

/** Get the Monte Carlo Barostat frequency in time speps */
int OpenMMFrEnergyST::getMCBarostatFrequency(void)
{
Expand Down
4 changes: 4 additions & 0 deletions corelib/src/libs/SireMove/openmmfrenergyst.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -122,6 +122,9 @@ namespace SireMove
bool getMCBarostat(void);
void setMCBarostat(bool);

bool getMCBarostatMembrane(void);
void setMCBarostatMembrane(bool);

void setMCBarostatFrequency(int);
int getMCBarostatFrequency(void);

Expand Down Expand Up @@ -247,6 +250,7 @@ namespace SireMove
double Andersen_frequency;

bool MCBarostat_flag;
bool MCBarostat_membrane_flag;
int MCBarostat_frequency;

QString ConstraintType;
Expand Down
Loading

0 comments on commit a5278e9

Please sign in to comment.