Skip to content

Commit

Permalink
Implement vdw1 and vdw2 lambdas
Browse files Browse the repository at this point in the history
  • Loading branch information
@martimunicoy
martimunicoy committed Feb 12, 2024
1 parent 692ebfd commit 36367de
Showing 1 changed file with 169 additions and 9 deletions.
178 changes: 169 additions & 9 deletions peleffy/topology/alchemistry.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -162,7 +162,8 @@ def connections(self):

def get_alchemical_topology(self, fep_lambda=None, coul_lambda=None,
coul1_lambda=None, coul2_lambda=None,
vdw_lambda=None, bonded_lambda=None):
vdw_lambda=None, vdw1_lambda=None,
vdw2_lambda=None, bonded_lambda=None):
"""
Given a lambda, it returns an alchemical topology after
combining both input topologies.
Expand Down Expand Up @@ -195,6 +196,20 @@ def get_alchemical_topology(self, fep_lambda=None, coul_lambda=None,
affects van der Waals parameters. It needs to be contained
between 0 and 1. It has precedence over fep_lambda.
Default is None
vdw1_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 1. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 1. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
vdw2_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 2. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 2. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
bonded_lambda : float
The value to define a coulombic lambda. This lambda only
affects bonded parameters. It needs to be contained
Expand All @@ -212,10 +227,13 @@ def get_alchemical_topology(self, fep_lambda=None, coul_lambda=None,
coul1_lambda = Coulombic1Lambda(coul1_lambda)
coul2_lambda = Coulombic2Lambda(coul2_lambda)
vdw_lambda = VanDerWaalsLambda(vdw_lambda)
vdw1_lambda = VanDerWaals1Lambda(vdw1_lambda)
vdw2_lambda = VanDerWaals2Lambda(vdw2_lambda)
bonded_lambda = BondedLambda(bonded_lambda)

lambda_set = LambdaSet(fep_lambda, coul_lambda, coul1_lambda,
coul2_lambda, vdw_lambda, bonded_lambda)
coul2_lambda, vdw_lambda, vdw1_lambda,
vdw2_lambda, bonded_lambda)

alchemical_topology = self.topology_from_lambda_set(lambda_set)

Expand Down Expand Up @@ -251,7 +269,7 @@ def topology_from_lambda_set(self, lambda_set):
atom.apply_lambda(["sigma", "epsilon", "born_radius",
"SASA_radius", "nonpolar_gamma",
"nonpolar_alpha"],
lambda_set.get_lambda_for_vdw(),
lambda_set.get_lambda_for_vdw1(),
reverse=False)
atom.apply_lambda(["charge"],
lambda_set.get_lambda_for_coulomb1(),
Expand All @@ -261,7 +279,7 @@ def topology_from_lambda_set(self, lambda_set):
atom.apply_lambda(["sigma", "epsilon", "born_radius",
"SASA_radius", "nonpolar_gamma",
"nonpolar_alpha"],
lambda_set.get_lambda_for_vdw(),
lambda_set.get_lambda_for_vdw2(),
reverse=True)
atom.apply_lambda(["charge"],
lambda_set.get_lambda_for_coulomb2(),
Expand Down Expand Up @@ -873,6 +891,7 @@ def molecule2_to_pdb(self, path):
def rotamer_library_to_file(self, path, fep_lambda=None,
coul_lambda=None, coul1_lambda=None,
coul2_lambda=None, vdw_lambda=None,
vdw1_lambda=None, vdw2_lambda=None,
bonded_lambda=None):
"""
It saves the alchemical rotamer library, which is the combination
Expand Down Expand Up @@ -909,6 +928,20 @@ def rotamer_library_to_file(self, path, fep_lambda=None,
affects van der Waals parameters. It needs to be contained
between 0 and 1. It has precedence over fep_lambda.
Default is None
vdw1_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 1. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 1. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
vdw2_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 2. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 2. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
bonded_lambda : float
The value to define a coulombic lambda. This lambda only
affects bonded parameters. It needs to be contained
Expand All @@ -927,14 +960,19 @@ def rotamer_library_to_file(self, path, fep_lambda=None,
coul1_lambda = Coulombic1Lambda(coul1_lambda)
coul2_lambda = Coulombic2Lambda(coul2_lambda)
vdw_lambda = VanDerWaalsLambda(vdw_lambda)
vdw1_lambda = VanDerWaals1Lambda(vdw1_lambda)
vdw2_lambda = VanDerWaals2Lambda(vdw2_lambda)
bonded_lambda = BondedLambda(bonded_lambda)

lambda_set = LambdaSet(fep_lambda, coul_lambda, coul1_lambda,
coul2_lambda, vdw_lambda, bonded_lambda)
coul2_lambda, vdw_lambda, vdw1_lambda,
vdw2_lambda, bonded_lambda)

if (at_least_one and
lambda_set.get_lambda_for_bonded() == 0.0 and
lambda_set.get_lambda_for_vdw() == 0.0 and
lambda_set.get_lambda_for_vdw1() == 0.0 and
lambda_set.get_lambda_for_vdw2() == 0.0 and
lambda_set.get_lambda_for_coulomb() == 0.0 and
lambda_set.get_lambda_for_coulomb1() == 0.0 and
lambda_set.get_lambda_for_coulomb2() == 0.0):
Expand All @@ -944,6 +982,8 @@ def rotamer_library_to_file(self, path, fep_lambda=None,
elif (at_least_one and
lambda_set.get_lambda_for_bonded() == 1.0 and
lambda_set.get_lambda_for_vdw() == 1.0 and
lambda_set.get_lambda_for_vdw1() == 1.0 and
lambda_set.get_lambda_for_vdw2() == 1.0 and
lambda_set.get_lambda_for_coulomb() == 1.0 and
lambda_set.get_lambda_for_coulomb1() == 1.0 and
lambda_set.get_lambda_for_coulomb2() == 1.0):
Expand Down Expand Up @@ -988,6 +1028,7 @@ def rotamer_library_to_file(self, path, fep_lambda=None,
def obc_parameters_to_file(self, path, fep_lambda=None,
coul_lambda=None, coul1_lambda=None,
coul2_lambda=None, vdw_lambda=None,
vdw1_lambda=None, vdw2_lambda=None,
bonded_lambda=None):
"""
It saves the alchemical OBC parameters, which is the combination
Expand Down Expand Up @@ -1030,6 +1071,20 @@ def obc_parameters_to_file(self, path, fep_lambda=None,
affects van der Waals parameters. It needs to be contained
between 0 and 1. It has precedence over fep_lambda.
Default is None
vdw1_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 1. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 1. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
vdw2_lambda : float
The value to define a vdw lambda for exclusive atoms
of molecule 2. This lambda only affects van der Waals
parameters of exclusive atoms of molecule 2. It needs
to be contained between 0 and 1.
It has precedence over vdw_lambda and fep_lambda.
Default is None
bonded_lambda : float
The value to define a coulombic lambda. This lambda only
affects bonded parameters. It needs to be contained
Expand All @@ -1055,10 +1110,13 @@ def obc_parameters_to_file(self, path, fep_lambda=None,
coul1_lambda = Coulombic1Lambda(coul1_lambda)
coul2_lambda = Coulombic2Lambda(coul2_lambda)
vdw_lambda = VanDerWaalsLambda(vdw_lambda)
vdw1_lambda = VanDerWaals1Lambda(vdw1_lambda)
vdw2_lambda = VanDerWaals2Lambda(vdw2_lambda)
bonded_lambda = BondedLambda(bonded_lambda)

lambda_set = LambdaSet(fep_lambda, coul_lambda, coul1_lambda,
coul2_lambda, vdw_lambda, bonded_lambda)
coul2_lambda, vdw_lambda, vdw1_lambda,
vdw2_lambda, bonded_lambda)

# Define mappers
mol1_mapped_atoms = [atom_pair[0] for atom_pair in self.mapping]
Expand Down Expand Up @@ -1257,6 +1315,20 @@ class VanDerWaalsLambda(Lambda):
_TYPE = "vdw"


class VanDerWaals1Lambda(Lambda):
"""
It defines the VanDerWaalsLambda1 class. It affects only van der Waals
parameters involving exclusive atoms of molecule 1.
"""
_TYPE = "vdw1"

class VanDerWaals2Lambda(Lambda):
"""
It defines the VanDerWaalsLambda1 class. It affects only van der Waals
parameters involving exclusive atoms of molecule 2.
"""
_TYPE = "vdw2"

class BondedLambda(Lambda):
"""
It defines the BondedLambda class. It affects only bonded parameters.
Expand All @@ -1270,7 +1342,7 @@ class LambdaSet(object):
"""

def __init__(self, fep_lambda, coul_lambda, coul1_lambda, coul2_lambda,
vdw_lambda, bonded_lambda):
vdw_lambda, vdw1_lambda, vdw2_lambda, bonded_lambda):
"""
It initializes a LambdaSet object which stores all the different
types of lambda.
Expand All @@ -1287,6 +1359,10 @@ def __init__(self, fep_lambda, coul_lambda, coul1_lambda, coul2_lambda,
The coulombic lambda for exclusive atoms of molecule 2
vdw_lambda : a peleffy.topology.alchemy.VanDerWaalsLambda object
The van der Waals lambda
vdw1_lambda : a peleffy.topology.alchemy.VanDerWaals1Lambda object
The van der Waals lambda for exclusive atoms of molecule 1
vdw2_lambda : a peleffy.topology.alchemy.VanDerWaals2Lambda object
The van der Waals lambda for exclusive atoms of molecule 2
bonded_lambda : a peleffy.topology.alchemy.BondedLambda object
The bonded lambda
"""
Expand All @@ -1308,6 +1384,12 @@ def __init__(self, fep_lambda, coul_lambda, coul1_lambda, coul2_lambda,
if not isinstance(vdw_lambda,
peleffy.topology.alchemistry.VanDerWaalsLambda):
raise TypeError('Invalid vdw_lambda supplied to LambdaSet')
if not isinstance(vdw1_lambda,
peleffy.topology.alchemistry.VanDerWaals1Lambda):
raise TypeError('Invalid vdw_lambda1 supplied to LambdaSet')
if not isinstance(vdw2_lambda,
peleffy.topology.alchemistry.VanDerWaals2Lambda):
raise TypeError('Invalid vdw2_lambda supplied to LambdaSet')
if not isinstance(bonded_lambda,
peleffy.topology.alchemistry.BondedLambda):
raise TypeError('Invalid bonded_lambda supplied to LambdaSet')
Expand All @@ -1317,6 +1399,8 @@ def __init__(self, fep_lambda, coul_lambda, coul1_lambda, coul2_lambda,
self._coul1_lambda = coul1_lambda
self._coul2_lambda = coul2_lambda
self._vdw_lambda = vdw_lambda
self._vdw1_lambda = vdw1_lambda
self._vdw2_lambda = vdw2_lambda
self._bonded_lambda = bonded_lambda

@property
Expand Down Expand Up @@ -1379,6 +1463,30 @@ def vdw_lambda(self):
"""
return self._vdw_lambda

@property
def vdw1_lambda(self):
"""
It returns the vdw1_lambda value.
Returns
-------
vdw1_lambda : float
The value of the vdw1_lambda
"""
return self._vdw1_lambda

@property
def vdw2_lambda(self):
"""
It returns the vdw2_lambda value.
Returns
-------
vdw2_lambda : float
The value of the vdw2_lambda
"""
return self._vdw2_lambda

@property
def bonded_lambda(self):
"""
Expand All @@ -1393,12 +1501,14 @@ def bonded_lambda(self):

def get_lambda_for_vdw(self):
"""
It returns the lambda to be applied on van der Waals parameters.
It returns the lambda to be applied on van der Waals parameters of
both molecules.
Returns
-------
lambda_value : float
The lambda value to be applied on van der Waals parameters
The lambda value to be applied on van der Waals parameters of
both molecules
"""
if self.vdw_lambda.is_set:
lambda_value = self.vdw_lambda.value
Expand All @@ -1411,6 +1521,56 @@ def get_lambda_for_vdw(self):

return lambda_value

def get_lambda_for_vdw1(self):
"""
It returns the lambda to be applied on van der Waals parameters of
exclusive atoms of molecule 1.
Returns
-------
lambda_value : float
The lambda value to be applied on van der Waals parameters of
exclusive atoms of molecule 1
"""
if self.vdw1_lambda.is_set:
lambda_value = self.vdw1_lambda.value

elif self.vdw_lambda.is_set:
lambda_value = self.vdw_lambda.value

elif self.fep_lambda.is_set:
lambda_value = self.fep_lambda.value

else:
lambda_value = 0.0

return lambda_value

def get_lambda_for_vdw2(self):
"""
It returns the lambda to be applied on van der Waals parameters of
exclusive atoms of molecule 2.
Returns
-------
lambda_value : float
The lambda value to be applied on van der Waals parameters of
exclusive atoms of molecule 2
"""
if self.vdw2_lambda.is_set:
lambda_value = self.vdw2_lambda.value

elif self.vdw_lambda.is_set:
lambda_value = self.vdw_lambda.value

elif self.fep_lambda.is_set:
lambda_value = self.fep_lambda.value

else:
lambda_value = 0.0

return lambda_value

def get_lambda_for_coulomb(self):
"""
It returns the lambda to be applied on Coulomb parameters of
Expand Down

0 comments on commit 36367de

Please sign in to comment.