Skip to content

Commit

Permalink
monte carlo prototyping and docs
Browse files Browse the repository at this point in the history
  • Loading branch information
mrbuche committed Aug 18, 2022
1 parent 19b6a2d commit 89156f3
Show file tree
Hide file tree
Showing 3 changed files with 307 additions and 52 deletions.
44 changes: 28 additions & 16 deletions statmechcrack/isometric.py
Original file line number Diff line number Diff line change
Expand Up @@ -130,13 +130,14 @@ def beta_A_abs_isometric(self, v, approach='asymptotic'):
numpy.ndarray: The absolute nondimensional Helmholtz free energy.
"""
s_hat = self.minimize_beta_U(v)[2]
lambda_hat = s_hat[-self.M:]
return 0.5*self.M*np.log(
self.alpha**2*self.varepsilon/np.pi
) + self.beta_A_0_isometric(v, lambda_hat)
if approach == 'asymptotic':
s_hat = self.minimize_beta_U(v)[2]
lambda_hat = s_hat[-self.M:]
return 0.5*self.M*np.log(
self.alpha**2*self.varepsilon/np.pi
) + self.beta_A_0_isometric(v, lambda_hat)

def beta_A_isometric(self, v, approach='asymptotic'):
def beta_A_isometric(self, v, approach='asymptotic', **kwargs):
r"""The relative nondimensional Helmholtz free energy
as a function of the nondimensional end separation,
Expand All @@ -147,6 +148,8 @@ def beta_A_isometric(self, v, approach='asymptotic'):
v (array_like): The nondimensional end separation.
approach (str, optional, default='asymptotic'):
The calculation approach.
**kwargs: Arbitrary keyword arguments.
Passed to :meth:`~.beta_A_isometric_monte_carlo`.
Returns:
numpy.ndarray: The relative nondimensional Helmholtz free energy.
Expand Down Expand Up @@ -180,8 +183,13 @@ def beta_A_isometric(self, v, approach='asymptotic'):
>>> plt.show()
"""
return self.beta_A_abs_isometric(v, approach=approach) \
- self.beta_A_abs_isometric(1, approach=approach)
if approach == 'asymptotic':
beta_A_isometric = \
self.beta_A_abs_isometric(v, approach=approach) - \
self.beta_A_abs_isometric(1, approach=approach)
elif approach == 'monte carlo':
beta_A_isometric = self.beta_A_isometric_monte_carlo(v, **kwargs)
return beta_A_isometric

def p_isometric(self, v, approach='asymptotic', **kwargs):
r"""The nondimensional end force
Expand Down Expand Up @@ -237,7 +245,7 @@ def p_isometric(self, v, approach='asymptotic', **kwargs):
p_isometric = self.p_isometric_monte_carlo(v, **kwargs)
return p_isometric

def k_isometric(self, v, approach='asymptotic'):
def k_isometric(self, v, approach='asymptotic', **kwargs):
r"""The nondimensional forward reaction rate coefficient
as a function of the nondimensional end separation
in the isometric ensemble.
Expand All @@ -246,6 +254,8 @@ def k_isometric(self, v, approach='asymptotic'):
v (array_like): The nondimensional end separation.
approach (str, optional, default='asymptotic'):
The calculation approach.
**kwargs: Arbitrary keyword arguments.
Passed to :meth:`~.k_isometric_monte_carlo`.
Returns:
numpy.ndarray: The nondimensional forward reaction rate.
Expand Down Expand Up @@ -277,13 +287,15 @@ def k_isometric(self, v, approach='asymptotic'):
>>> plt.show()
"""
return (
self.Q_isometric(v, approach=approach, transition_state=True) /
self.Q_isometric(v, approach=approach)
) / (
self.Q_isometric(1, approach=approach, transition_state=True) /
self.Q_isometric(1, approach=approach)
)
if approach == 'asymptotic':
k_isometric = (
self.Q_isometric(v, transition_state=True)/self.Q_isometric(v)
) / (
self.Q_isometric(1, transition_state=True)/self.Q_isometric(1)
)
elif approach == 'monte carlo':
k_isometric = self.k_isometric_monte_carlo(v, **kwargs)
return k_isometric

def Q_0_isometric(self, v, lambda_):
r"""The nondimensional isometric partition function
Expand Down
40 changes: 26 additions & 14 deletions statmechcrack/isotensional.py
Original file line number Diff line number Diff line change
Expand Up @@ -125,11 +125,12 @@ def beta_G_abs_isotensional(self, p, approach='asymptotic'):
numpy.ndarray: The absolute nondimensional Gibbs free energy.
"""
s_hat = self.minimize_beta_Pi(p)[2]
lambda_hat = s_hat[-self.M:]
return self.beta_G_0_abs_isotensional(p, lambda_hat)
if approach == 'asymptotic':
s_hat = self.minimize_beta_Pi(p)[2]
lambda_hat = s_hat[-self.M:]
return self.beta_G_0_abs_isotensional(p, lambda_hat)

def beta_G_isotensional(self, p, approach='asymptotic'):
def beta_G_isotensional(self, p, approach='asymptotic', **kwargs):
r"""The relative nondimensional Gibbs free energy
as a function of the nondimensional end force,
Expand All @@ -140,6 +141,8 @@ def beta_G_isotensional(self, p, approach='asymptotic'):
p (array_like): The nondimensional end force.
approach (str, optional, default='asymptotic'):
The calculation approach.
**kwargs: Arbitrary keyword arguments.
Passed to :meth:`~.beta_G_isotensional_monte_carlo`.
Returns:
numpy.ndarray: The relative nondimensional Gibbs free energy.
Expand Down Expand Up @@ -174,8 +177,13 @@ def beta_G_isotensional(self, p, approach='asymptotic'):
>>> plt.show()
"""
return self.beta_G_abs_isotensional(p, approach=approach) \
- self.beta_G_abs_isotensional(0, approach=approach)
if approach == 'asymptotic':
G_isotensional = \
self.beta_G_abs_isotensional(p, approach=approach) - \
self.beta_G_abs_isotensional(0, approach=approach)
elif approach == 'monte carlo':
G_isotensional = self.beta_G_isotensional_monte_carlo(p, **kwargs)
return G_isotensional

def v_isotensional(self, p, approach='asymptotic', **kwargs):
r"""The nondimensional end separation
Expand Down Expand Up @@ -237,7 +245,7 @@ def v_isotensional(self, p, approach='asymptotic', **kwargs):
elif approach == 'monte carlo':
return self.v_isotensional_monte_carlo(p, **kwargs)

def k_isotensional(self, p, approach='asymptotic'):
def k_isotensional(self, p, approach='asymptotic', **kwargs):
r"""The nondimensional forward reaction rate coefficient
as a function of the nondimensional end force
in the isotensional ensemble.
Expand Down Expand Up @@ -278,13 +286,17 @@ def k_isotensional(self, p, approach='asymptotic'):
>>> plt.show()
"""
return (
self.Z_isotensional(p, approach=approach, transition_state=True) /
self.Z_isotensional(p, approach=approach)
) / (
self.Z_isotensional(0, approach=approach, transition_state=True) /
self.Z_isotensional(0, approach=approach)
)
if approach == 'asymptotic':
k_isotensional = (
self.Z_isotensional(p, transition_state=True) /
self.Z_isotensional(p)
) / (
self.Z_isotensional(0, transition_state=True) /
self.Z_isotensional(0)
)
elif approach == 'monte carlo':
k_isotensional = self.k_isotensional_monte_carlo(p, **kwargs)
return k_isotensional

def Z_0_isotensional(self, p, lambda_):
r"""The nondimensional isotensional partition function
Expand Down
Loading

0 comments on commit 89156f3

Please sign in to comment.