"""A module for handling potentials.This module contains several different classes representing potentials,each having methods to compute relevant nondimensional quantities asfunctions of nondimensional force or stretch.This module also contains the parent class ``Potential`` that is used toassign a potential to a given model using keyword arguments.Examples: Create a Lennard-Jones potential model with a nondimensional potential energy scale of 8 and evaluate the nondimensional potential energy at a stretch of 1.23: >>> from ufjc.potential import LennardJonesPotential >>> model = LennardJonesPotential(varepsilon=8) >>> model.beta_u(1.23) 4.046654314368616 Do the same with the Lenard-Jones-FENE potential: >>> from ufjc.potential import LJFENEPotential >>> model = LJFENEPotential(varepsilon=(8, 8)) >>> model.beta_u(1.23) 8.510502022381505 Create a single-link model in one dimension, instantiate it with the Morse potential, and compute the incremental link stretch under a nondimensional force of 8: >>> from ufjc.potential import Potential >>> class Link1D(Potential): ... def __init__(self, **kwargs): ... Potential.__init__(self, **kwargs) >>> _ = Link1D(potential='morse').delta_lambda(8) >>> Link1D(potential='lj-fene').eta_link(1) 184.0"""# Import external modulesimportnumpyasnpfromscipy.specialimportlambertw
[docs]classHarmonicPotential(object):r"""The harmonic potential. Attributes: varepsilon (float): The nondimensional energy scale. kappa (float): The nondimensional stiffness :math:`\kappa=\varepsilon`. c (float): The correction parameter :math:`c=1`. """def__init__(self,**kwargs):"""Initializes the ``HarmonicPotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88). """self.varepsilon=kwargs.get('varepsilon',88)self.kappa=self.varepsilonself.c=1
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \frac{1}{2}(\lambda-1)^2 . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """return0.5*(lambda_-1)**2
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = \varepsilon(\lambda - 1) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). Example: Compute the nondimensional force at a sample stretch: >>> from ufjc.potential import HarmonicPotential >>> HarmonicPotential().eta_link(1.8) 70.4 """returnself.varepsilon*(lambda_-1)
[docs]defdelta_lambda(self,eta):r"""The incremental stretch as a function of nondimensional force, .. math:: \Delta\lambda(\eta) = \frac{\eta}{\varepsilon} . Args: eta (array_like): The nondimensional force(s). Returns: numpy.ndarray: The incremental stretch(s). """returneta/self.varepsilon
[docs]classLogSquaredPotential(object):r"""The log-squared potential :cite:`potential-mao2017rupture`. Attributes: varepsilon (float): The nondimensional energy scale. kappa (float): The nondimensional stiffness :math:`\kappa=\varepsilon`. c (float): The correction parameter :math:`c=2/5`. eta_max (float): The maximum nondimensional force :math:`\eta_\mathrm{max} = e^{-1}\varepsilon`. lambda_max (float): The stretch at the maximum nondimensional force, :math:`\lambda_\mathrm{max} = e^{1}`. """def__init__(self,**kwargs):"""Initializes the ``LogSquaredPotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88). """self.varepsilon=kwargs.get('varepsilon',88)self.kappa=self.varepsilonself.c=2/5self.eta_max=self.varepsilon/np.exp(1)self.lambda_max=np.exp(1)
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \frac{1}{2}\big[\ln(\lambda)\big]^2 . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """return0.5*np.log(lambda_)**2
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = \varepsilon\,\frac{\ln(\lambda)}{\lambda} . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). Example: Compute the nondimensional force at a sample stretch: >>> from ufjc.potential import LogSquaredPotential >>> LogSquaredPotential().eta_link(1.8) 28.736236950770266 """returnself.varepsilon*np.log(lambda_)/lambda_
[docs]defdelta_lambda(self,eta):r"""The incremental stretch as a function of nondimensional force, .. math:: \Delta\lambda(\eta) = e^{-\mathcal{W}(-\eta/\varepsilon)} ,\qquad \eta\in[0,\eta_\mathrm{max}] . Args: eta (array_like): The nondimensional force(s). Returns: numpy.ndarray: The incremental stretch(s). """return(np.exp(-lambertw(-eta/self.varepsilon))-1).real
[docs]classMorsePotential(object):r"""The Morse potential :cite:`potential-morse1929diatomic`. Attributes: varepsilon (float): The nondimensional energy scale. alpha (float): The Morse parameter. kappa (float): The nondimensional stiffness :math:`\kappa=2\varepsilon\alpha^2`. c (float): The correction parameter :math:`c=1/(1+3\alpha/2)`. eta_max (float): The maximum nondimensional force :math:`\eta_\mathrm{max} = \sqrt{\kappa\varepsilon/8}`. lambda_max (float): The stretch at the maximum nondimensional force, :math:`\lambda_\mathrm{max} = 1+\ln(2)/\alpha`. """def__init__(self,**kwargs):"""Initializes the ``MorsePotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88) and ``alpha`` (default 1). """self.varepsilon=kwargs.get('varepsilon',88)self.alpha=kwargs.get('alpha',1)self.kappa=2*self.varepsilon*self.alpha**2self.c=1/(1+3/2*self.alpha)self.eta_max=np.sqrt(self.kappa*self.varepsilon/8)self.lambda_max=1+np.log(2)/self.alpha
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \left[1 - e^{-\alpha(\lambda - 1)}\right]^2 . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """return(1-np.exp(-self.alpha*(lambda_-1)))**2
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = 2\alpha\varepsilon e^{-\alpha(\lambda - 1)} \left[1 - e^{-\alpha(\lambda - 1)}\right]^2 . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). Example: Compute the nondimensional force at a sample stretch: >>> from ufjc.potential import MorsePotential >>> _ = MorsePotential().eta_link(1.23) """return2*self.alpha*self.varepsilon* \
np.exp(-self.alpha*(lambda_-1))* \
(1-np.exp(-self.alpha*(lambda_-1)))
[docs]defdelta_lambda(self,eta):r"""The incremental stretch as a function of nondimensional force, .. math:: \Delta\lambda(\eta) = \ln\left(\frac{2}{1 + \sqrt{1 - \eta/\eta_\mathrm{max}}}\right)^{1/\alpha} ,\qquad \eta\in[0,\eta_\mathrm{max}] . Args: eta (array_like): The nondimensional force(s). Returns: numpy.ndarray: The incremental stretch(s). """returnnp.log(2/(1+np.sqrt(1-eta/self.eta_max)))/self.alpha
[docs]classLennardJonesPotential(object):r"""The Lennard-Jones potential :cite:`potential-jones1924determinationii`. Attributes: varepsilon (float): The nondimensional energy scale. kappa (float): The nondimensional stiffness :math:`\kappa=72\varepsilon`. c (float): The correction parameter :math:`c=2/23`. eta_max (float): The maximum nondimensional force :math:`\eta_\mathrm{max} = \eta(\lambda_\mathrm{max})`. lambda_max (float): The stretch at the maximum nondimensional force, :math:`\lambda_\mathrm{max} = (13/7)^{1/6}`. """def__init__(self,**kwargs):"""Initializes the ``LennardJonesPotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88). """self.varepsilon=kwargs.get('varepsilon',88)self.kappa=72*self.varepsilonself.c=2/23self.lambda_max=(13/7)**(1/6)self.eta_max=self.eta_link(self.lambda_max)
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \frac{1}{\lambda^{12}} - \frac{2}{\lambda^6} + 1 . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """return1/lambda_**12-2/lambda_**6+1
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = 12\varepsilon\left( \frac{1}{\lambda^7} - \frac{1}{\lambda^{13}}\right) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). """returnself.varepsilon*(12/lambda_**7-12/lambda_**13)
[docs]classLJFENEPotential(object):r"""The Lennard-Jones-FENE potential :cite:`potential-kremer1990dynamics`. Attributes: varepsilon (float): The nondimensional energy scale. kappa (float): The nondimensional stiffness. c (float): The correction parameter. eta_max (float): The maximum nondimensional force :math:`\eta_\mathrm{max} = \eta(\lambda_\mathrm{max})`. lambda_max (float): The stretch at the maximum nondimensional force. """def__init__(self,**kwargs):"""Initializes the ``LJFENEPotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default (88, 230)) and ``lambda_max`` (default 1.5). """self.varepsilon_1,self.varepsilon_2= \
kwargs.get('varepsilon',(88,230))self.varepsilon=self.varepsilon_2self.lambda_max=kwargs.get('lambda_max',1.5)self.kappa=72*self.varepsilon_1+self.varepsilon_2* \
(self.lambda_max**2+1)/(self.lambda_max**2-1)**2self.c=1/(1-(-1512+(6*self.lambda_max**2+2)/(self.lambda_max**2-1)**3)/(2*self.kappa/self.varepsilon_2))self.eta_max=self.eta_link(self.lambda_max)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = 12\varepsilon_1\left( \frac{1}{\lambda^7} - \frac{1}{\lambda^{13}} \right) + \frac{\varepsilon_2\lambda} {\lambda_\mathrm{max}^2 - \lambda^2} . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). """lambda_fene=(lambda_<self.lambda_max)*lambda_returnself.varepsilon_1*(12/lambda_**7-12/lambda_**13)+ \
self.varepsilon_2*lambda_fene/(self.lambda_max**2-lambda_fene**2)*(lambda_<self.lambda_max)
[docs]classMiePotential(object):r"""The Mie potential :cite:`potential-mie1903kinetischen`. Attributes: varepsilon (float): The nondimensional energy scale. n (float): The repulsive exponent. m (float): The attractive exponent. kappa (float): The nondimensional stiffness :math:`\kappa=nm\varepsilon`. c (float): The correction parameter :math:`c=\frac{4m(m+1)-2n(n+1)}{2m(m^2+5m+4)-n(n^2+5n+4)}`. eta_max (float): The maximum nondimensional force :math:`\eta_\mathrm{max} = \eta(\lambda_\mathrm{max})`. lambda_max (float): The stretch at the maximum nondimensional force, :math:`\lambda_\mathrm{max}=\left[\frac{n+1}{m+1}\right]^{1/(n-m)}` . """def__init__(self,**kwargs):"""Initializes the ``MiePotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88) ``n`` (default 12) and ``m`` (default 6). """self.varepsilon=kwargs.get('varepsilon',88)self.m=kwargs.get('m',6)self.n=kwargs.get('n',12)self.kappa=self.n*self.m*self.varepsilonself.c=(4*self.m*(self.m+1)-2*self.n*(self.n+1))/(2*self.m*(self.m**2+5*self.m+4)-self.n*(self.n**2+5*self.n+4))self.lambda_max=((self.n+1)/(self.m+1))**(1/(self.n-self.m))self.eta_max=self.eta_link(self.lambda_max)
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \frac{1}{(n-m)} \left(\frac{m}{\lambda^n} - \frac{n}{\lambda^m}\right) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """return(self.m/lambda_**self.n-self.n/lambda_**self.m)/(self.n-self.m)
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = \varepsilon\, \frac{nm}{(n-m)} \left(\frac{1}{\lambda^{m+1}} - \frac{1}{\lambda^{n+1}}\right) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). """returnself.varepsilon*self.n*self.m/(self.n-self.m)*(1/lambda_**(self.m+1)-1/lambda_**(self.n+1))
[docs]classPolynomialPotential(object):r"""A polynomial potential. Attributes: varepsilon (float): The nondimensional energy scale. coefficients (array_like): The coefficients :math:`a_k`. kappa (float): The nondimensional stiffness :math:`\kappa=a_1\varepsilon`. c (float): The correction parameter :math:`c=(1 - \frac{a_2}{2a_1})^{-1}`. """def__init__(self,**kwargs):"""Initializes the ``PolynomialPotential`` class. Args: **kwargs: Arbitrary keyword arguments. Can be used to specify ``varepsilon`` (default 88) and ``coefficients`` (default 1). """self.varepsilon=kwargs.get('varepsilon',88)coef=np.array(kwargs.get('coefficients',[1,0]))self.eta_c=np.append(np.array([0]),coef)self.phi_c=np.append([0],self.eta_c)*np.append([0,0],[1/nforninrange(2,len(coef)+2)])self.kappa=self.varepsilon*self.eta_c[1]iflen(self.phi_c)>3:self.c=1/(1-self.phi_c[3]/self.phi_c[2]/2)
[docs]defphi(self,lambda_):r"""The scaled nondimensional potential energy function, .. math:: \phi(\lambda) = \sum_{k=2} \frac{a_{k-1}}{k}\, \lambda^k . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The scaled nondimensional potential energy(s). """returnnp.polynomial.Polynomial(self.phi_c)(lambda_-1)
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]defeta_link(self,lambda_):r"""The nondimensional force as a function of stretch, .. math:: \eta(\lambda) = \varepsilon\sum_{k=1} a_k\, \lambda^k . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional force(s). """returnself.varepsilon* \
np.polynomial.Polynomial(self.eta_c)(lambda_-1)
[docs]classCustomPotential(object):r"""A custom user-defined potential. Attributes: potential (str): The potential name. varepsilon (float): The nondimensional energy scale. phi (function): The scaled nondimensional energy function. eta_link (function): The nondimensional force as a function of stretch. delta_lambda (function): The incremental stretch as a function of the nondimensional force (optional). kappa (float): The nondimensional stiffness. c (float): The correction parameter. """def__init__(self,**kwargs):"""Initializes the ``CustomPotential`` class. Args: **kwargs: Arbitrary keyword arguments. """self.potential=kwargs.get('potential','custom')self.varepsilon=kwargs.get('varepsilon')self.phi=kwargs.get('phi')self.eta_link=kwargs.get('eta_link')self.kappa=kwargs.get('kappa')self.c=kwargs.get('c')# No defaultsif'delta_lambda'inkwargs:self.delta_lambda=kwargs.get('delta_lambda')if'lambda_max'inkwargs:self.lambda_max=kwargs.get('lambda_max')if'eta_max'inkwargs:self.eta_max=kwargs.get('eta_max')
[docs]defbeta_u(self,lambda_):r"""The nondimensional potential energy function, .. math:: \beta u(\lambda) = \varepsilon\phi(\lambda) . Args: lambda_ (array_like): The stretch(s). Returns: numpy.ndarray: The nondimensional potential energy(s). """returnself.varepsilon*self.phi(lambda_)
[docs]classPotential(object):r"""A class to assign a potential to a given model through inheritance. Attributes: potential (str): The potential type. pot (object): The potential model instance. """def__init__(self,**kwargs):"""Initializes the ``Potential`` class. Args: **kwargs: Arbitrary keyword arguments. Used to specify the potential and then passed to an instantiation of that potential. Note: An improperly-specified potential will default to harmonic: >>> from ufjc import uFJC >>> model = uFJC(potential='blah') Potential "blah" invalid, defaulting to "harmonic" potential. """self.potential=kwargs.get('potential','harmonic')ifself.potential=='harmonic':self.pot=HarmonicPotential(**kwargs)elifself.potential=='log-squared':self.pot=LogSquaredPotential(**kwargs)elifself.potential=='morse':self.pot=MorsePotential(**kwargs)elifself.potential=='lennard-jones':self.pot=LennardJonesPotential(**kwargs)elifself.potential=='lj-fene':self.pot=LJFENEPotential(**kwargs)elifself.potential=='mie':self.pot=MiePotential(**kwargs)elifself.potential=='polynomial':self.pot=PolynomialPotential(**kwargs)elifself.potential=='custom':self.pot=CustomPotential(**kwargs)else:print('Potential "'+self.potential+'" invalid, defaulting to "harmonic" potential.')self.pot=HarmonicPotential(**kwargs)def__getattr__(self,attr):"""Inherit attributes and methods from the chosen potential. Note: When accessing an attribute of the model to which the potential is assigned, but no attribute is found, and ``AttributeError`` will be raised with respect to the ``pot`` instance, not the overall model instance. This is sub-optimal but typically acceptable, but still should be repaired in the future. """returngetattr(self.pot,attr)