Merge pull request #1 from sandialabs/dev

more efficient and accurate monte carlo approach

Author: mrbuche

statmechcrack/isometric.py

@@ -83,7 +83,7 @@ def Q_isometric(self, v, transition_state=False):
                 >>> from statmechcrack import CrackIsometric
                 >>> v = np.linspace(1, 11, 250)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [10, 25, 100, 800]:
+                >>> for varepsilon in [10, 25, 100, 1000]:
                 ...     model = CrackIsometric(varepsilon=varepsilon)
                 ...     r_Q = (model.Q_isometric(v)/model.Q_isometric(1)
                 ...     )**(model.N**3/3/model.kappa)
@@ -170,7 +170,7 @@ def beta_A_isometric(self, v, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsometric
                 >>> v = np.linspace(1, 11, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [10, 25, 100, 800]:
+                >>> for varepsilon in [10, 25, 100, 1000]:
                 ...     model = CrackIsometric(varepsilon=varepsilon)
                 ...     beta_A = \
                 ...         model.beta_A_isometric(v, approach='asymptotic')
@@ -225,7 +225,7 @@ def p_isometric(self, v, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsometric
                 >>> v = np.linspace(1, 11, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [10, 25, 100, 800]:
+                >>> for varepsilon in [10, 25, 100, 1000]:
                 ...     model = CrackIsometric(varepsilon=varepsilon)
                 ...     p = model.p_isometric(v, approach='asymptotic')
                 ...     _ = plt.plot(v - 1, model.N**3/3/model.kappa*p,
@@ -276,7 +276,7 @@ def k_isometric(self, v, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsometric
                 >>> v = np.linspace(1, 11, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [10, 25, 100, 800]:
+                >>> for varepsilon in [10, 25, 100, 1000]:
                 ...     model = CrackIsometric(varepsilon=varepsilon)
                 ...     _ = plt.semilogy(
                 ...         v - 1, model.k_isometric(v),
@@ -415,7 +415,7 @@ def beta_A_0_isometric(self, v, lambda_):
 
         """
         return self.beta_A_0_abs_isometric(v, lambda_) \
-            - self.beta_A_0_abs_isometric(lambda_[0], lambda_)
+            - self.beta_A_0_abs_isometric(1, lambda_)
 
     def p_0_isometric(self, v, lambda_):
         r"""The nondimensional end force
@@ -579,7 +579,7 @@ def beta_A_b_isometric(self, v, lambda_):
 
         """
         return self.beta_A_b_isometric_abs(v, lambda_) \
-            - self.beta_A_b_isometric_abs(lambda_[0], lambda_)
+            - self.beta_A_b_isometric_abs(1, lambda_)
 
     def p_b_isometric(self, v, lambda_):
         r"""The nondimensional end force

statmechcrack/isotensional.py

@@ -75,7 +75,7 @@ def Z_isotensional(self, p, transition_state=False):
                 >>> from statmechcrack import CrackIsotensional
                 >>> rp = np.linspace(0, 10, 250)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [50, 100, 800]:
+                >>> for varepsilon in [50, 100, 1000]:
                 ...     model = CrackIsotensional(varepsilon=varepsilon)
                 ...     p = 3*model.kappa/model.N**3*rp
                 ...     r_Z = (model.Z_isotensional(0)/model.Z_isotensional(p)
@@ -164,7 +164,7 @@ def beta_G_isotensional(self, p, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsotensional
                 >>> rp = np.linspace(0, 10, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [50, 100, 800]:
+                >>> for varepsilon in [50, 100, 1000]:
                 ...     model = CrackIsotensional(varepsilon=varepsilon)
                 ...     p = 3*model.kappa/model.N**3*rp
                 ...     beta_G = model.beta_G_isotensional(
@@ -220,7 +220,7 @@ def v_isotensional(self, p, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsotensional
                 >>> rp = np.linspace(0, 10, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [50, 100, 800]:
+                >>> for varepsilon in [10, 25, 100, 1000]:
                 ...     model = CrackIsotensional(varepsilon=varepsilon)
                 ...     p = 3*model.kappa/model.N**3*rp
                 ...     v = model.v_isotensional(p, approach='asymptotic')
@@ -268,7 +268,7 @@ def k_isotensional(self, p, approach='asymptotic', **kwargs):
                 >>> from statmechcrack import CrackIsotensional
                 >>> rp = np.linspace(0, 10, 33)
                 >>> _ = plt.figure()
-                >>> for varepsilon in [50, 100, 800]:
+                >>> for varepsilon in [50, 100, 1000]:
                 ...     model = CrackIsotensional(varepsilon=varepsilon)
                 ...     p = 3*model.kappa/model.N**3*rp
                 ...     _ = plt.semilogy(

statmechcrack/mechanical.py

@@ -147,7 +147,9 @@ def beta_u(self, lambda_):
             numpy.ndarray: The nondimensional potential energy.
 
         """
-        return self.varepsilon*(1 - np.exp(-self.alpha*(lambda_ - 1)))**2
+        return self.varepsilon*(
+            1 - np.exp(-self.alpha*(lambda_ - 1))
+        )**2 + 1e88*(lambda_ > self.lambda_TS)
 
     def beta_u_p(self, lambda_):
         r"""The first derivative of the potential energy of a single bond