statmechcrack.monte_carlo

A module for the Metropolis-Hastings Markov chain Monte Carlo method.

This module consist of the class CrackMonteCarlo which contains methods to compute crack model quantities using the Metropolis-Hastings Markov chain Monte Carlo method [1].

class CrackMonteCarlo(**kwargs)

Bases: CrackMechanical

The Metropolis-Hastings Markov chain Monte Carlo class.

beta_A_isometric_monte_carlo(v, **kwargs)

The relative nondimensional Helmholtz free energy as a function of the nondimensional end separation in the isometric ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[\beta\Delta A(v) = -\ln\left\langle e^{-\beta\Delta A_0} \right\rangle_\star,\]

where the relative free energy of the reference system is

\[\Delta A_0\equiv A_0(v,\boldsymbol{\lambda})-A_0(1,\boldsymbol{\lambda}),\]

and where the isometric ensemble average here is

\[\left\langle \phi \right\rangle_\star \equiv \frac{ \int d\lambda \ e^{-\beta A_\star(\boldsymbol{\lambda})} \, \phi(\boldsymbol{\lambda}) }{ \int d\lambda \ e^{-\beta A_\star(\boldsymbol{\lambda})} } ,\]

where the free energy measure here is

\[A_\star(\boldsymbol{\lambda})\equiv A_0(1,\boldsymbol{\lambda}) + U_1(\boldsymbol{\lambda}).\]

Parameters:

• v (array_like) – The nondimensional end separation.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The relative nondimensional Helmholtz free energy.

Return type:

numpy.ndarray

beta_A_isometric_monte_carlo_serial(v, init_config, num_samples=1000000, **kwargs)

Serial calculation for beta_A_isometric_monte_carlo().

Parameters:

• v (array_like) – The nondimensional end separation.

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The relative nondimensional Helmholtz free energy.

Return type:

numpy.ndarray

beta_G_isotensional_monte_carlo(p, **kwargs)

The relative nondimensional Gibbs free energy as a function of the nondimensional end force in the isometric ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[\beta\Delta G(v) = -\ln\left\langle e^{-\beta\Delta G_0} \right\rangle_\star,\]

where the relative free energy of the reference system is

\[\Delta G_0\equiv G_0(p,\boldsymbol{\lambda})-G_0(0,\boldsymbol{\lambda}),\]

and where the isotensional ensemble average here is

\[\left\langle \phi \right\rangle_\star \equiv \frac{ \int d\lambda \ e^{-\beta G_\star(\boldsymbol{\lambda})} \, \phi(\boldsymbol{\lambda}) }{ \int d\lambda \ e^{-\beta G_\star(\boldsymbol{\lambda})} } ,\]

where the free energy measure here is

\[G_\star(\boldsymbol{\lambda})\equiv G_0(0,\boldsymbol{\lambda}) + U_1(\boldsymbol{\lambda}).\]

Parameters:

• p (array_like) – The nondimensional end force.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The relative nondimensional Gibbs free energy.

Return type:

numpy.ndarray

beta_G_isotensional_monte_carlo_serial(p, init_config, num_samples=1000000, **kwargs)

Serial calculation for beta_G_isotensional_monte_carlo().

Parameters:

• p (array_like) – The nondimensional end force.

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The relative nondimensional Gibbs free energy.

Return type:

numpy.ndarray

burn_in(init_config, num_burns=10000, tol=inf, **kwargs)

Burn-in (throw away samples) until in the effective sampling region.

This method runs a Monte Carlo calculation without utilizing samples, instead throwing them away until a desired number have been thrown away or a desired tolerance has been reached. This is typically done to obtain a burned-in configuration to use as the initial configuration for an actual Monte Carlo calculation.

Parameters:

• init_config (numpy.ndarray) – The initial configuration.

• num_burns (int, optional, default=10000) – The number of samples to burn.

• tol (float, optional, default=np.inf) – The burn-in tolerance.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config_append().

Returns:

The final, burned-in configuration.

Return type:

numpy.ndarray

k_isometric_monte_carlo(v, **kwargs)

The nondimensional forward reaction rate coefficient as a function of the nondimensional end force in the isometric ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[\frac{k(v)}{k(1)} = \frac{\left\langle e^{-\beta\Delta A_0^\ddagger} \right\rangle_\star^\ddagger}{\left\langle e^{-\beta\Delta A_0} \right\rangle_\star}.\]

Parameters:

• v (array_like) – The nondimensional end separation.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The nondimensional forward reaction rate.

Return type:

numpy.ndarray

k_isometric_monte_carlo_serial(v, init_config, num_samples=1000000, **kwargs)

Serial calculation for k_isometric_monte_carlo().

Parameters:

• v (array_like) – The nondimensional end separation.

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The nondimensional forward reaction rate.

Return type:

numpy.ndarray

k_isotensional_monte_carlo(p, **kwargs)

The nondimensional forward reaction rate coefficient as a function of the nondimensional end force in the isotensional ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[\frac{k(p)}{k(0)} = \frac{\left\langle e^{-\beta\Delta G_0^\ddagger} \right\rangle_\star^\ddagger}{\left\langle e^{-\beta\Delta G_0} \right\rangle_\star}.\]

Parameters:

• p (array_like) – The nondimensional end force.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The nondimensional forward reaction rate.

Return type:

numpy.ndarray

k_isotensional_monte_carlo_serial(p, init_config, num_samples=1000000, **kwargs)

Serial calculation for k_isotensional_monte_carlo().

Parameters:

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The nondimensional forward reaction rate.

Return type:

numpy.ndarray

mh_next_config(prev_config, urng=<built-in method uniform of numpy.random.mtrand.RandomState object>, **kwargs)

Generates accepted configurations given a previous configuration.

This function returns an accepted configuration given a previous configuration. Trial configurations are generated and repeatedly rejected until accepted as the next configuration based on the Metropolis-Hastings algorithm.

Parameters:

• prev_config (numpy.ndarray) – The previous configuration.

• urng (object, optional, default=np.random.uniform) – The uniform random number generator.

• **kwargs – Arbitrary keyword arguments. Passed to trial_config().

Returns:

The accepted (next) configuration.

Return type:

numpy.ndarray

mh_next_config_append(prev_configs, **kwargs)

Append next configuration onto a history of previous configurations.

This method calls mh_next_config() and appends the resulting next configuration onto a given history of previous configurations.

Parameters:

• prev_configs (numpy.ndarray) – The previous history of configurations.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The updated history of configurations.

Return type:

numpy.ndarray

p_isometric_monte_carlo(v, **kwargs)

The nondimensional end force as a function of the nondimensional end separation in the isometric ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[p(v) = e^{\beta\Delta A} \left\langle p_0 \, e^{-\beta\Delta A_0} \right\rangle_\star.\]

Parameters:

• v (array_like) – The nondimensional end separation.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The nondimensional end force.

Return type:

float

p_isometric_monte_carlo_serial(v, init_config, num_samples=1000000, **kwargs)

Serial calculation for p_isometric_monte_carlo().

Parameters:

• v (array_like) – The nondimensional end separation.

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The nondimensional end force.

Return type:

float

parallel_calculation(serial_fun, init_config, **kwargs)

Monte Carlo calculation averaged over several parallel processes. This method performs several Monte Carlo calculations in parallel and returns the average result from all processes. Each Monte Carlo calculation is performed using serial_fun(). The default is to utilize all available processors.

Parameters:

• serial_fun (function) – The function for a single Monte Carlo calculation.

• init_config (np.ndarray) – The initial configuration for the burn-in process.

• **kwargs – Arbitrary keyword arguments. Passed to burn_in() and serial_fun().

Returns:

The process-averaged result.

Return type:

float

Note

To avoid issues with multiprocessing and MacOS/Windows, num_processes is set to 1 for anything but Linux.

trial_config(prev_config, cov_config=0.01)

Generates trial configurations given a previous configuration.

Parameters:

• prev_config (numpy.ndarray) – The previous configuration.

• cov_config (float, optional, default=1e-2) – The covariance.

• **kwargs – Arbitrary keyword arguments.

Returns:

The generated trial configuration.

Return type:

numpy.ndarray

v_isotensional_monte_carlo(p, **kwargs)

The nondimensional end separation as a function of the nondimensional end force in the isotensional ensemble, using a Metropolis-Hastings Markov chain Monte Carlo calculation of the ensemble average

\[v(p) = e^{\beta\Delta G} \left\langle v_0 \, e^{-\beta\Delta G_0} \right\rangle_\star.\]

Parameters:

• p (array_like) – The nondimensional end force.

• **kwargs – Arbitrary keyword arguments. Passed to parallel_calculation().

Returns:

The nondimensional end separation.

Return type:

float

v_isotensional_monte_carlo_serial(p, init_config, num_samples=1000000, **kwargs)

Serial calculation for v_isotensional_monte_carlo().

Parameters:

• init_config (np.ndarray) – The initial configuration, typically already burned-in.

• num_samples (int, optional, default=1000000) – The number of samples to use.

• **kwargs – Arbitrary keyword arguments. Passed to mh_next_config().

Returns:

The nondimensional end separation.

Return type:

float

---

References

[1]

James M Haile. Molecular Dynamics Simulation: Elementary Methods. Wiley, 1992. URL: https://www.wiley.com/Molecular+Dynamics+Simulation:+Elementary+Methods-p-9780471184393.