Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add support for OpenMM MonteCarloMembraneBarostat #281

Merged
merged 1 commit into from
Feb 4, 2025
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
45 changes: 37 additions & 8 deletions corelib/src/libs/SireMove/openmmfrenergydt.cpp
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
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
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
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
Loading